Pesticide Heat Mapping 7/29/24
DaSEH Project
This was started at the DaSEH Code-A-Thon! Should give some credit to them for coding help.
Analysis Plan from Thesis Proposal
Aim 3 will explore and describe these co-exposures to visually display the county-level spatial distribution of this co-exposure across the United States. Priority areas for public health interventions are areas that have high co-exposures and high vulnerability to environmental hazards like heat. The primary outcome variable for mapping is the co-exposure of heat and glyphosate, defined here as the heat-glyphosate index calculated as 85th percentile high temperature*glyphosate application per square kilometer. I will map this primary outcome on a scale from white to bright red where red is the highest co-exposure. I hypothesize that the Midwest and Great Plains regions will have the highest levels of co-exposure due to the large amounts of glyphosate resistant crops grown there.
Relevant exposures of public health and policy interest include the CDC’s Social Vulnerability Index (SVI), most prevalent type of glyphosate resistant crop, and urban pesticide usage. I will overlay SVI over the glyphosate-heat co-exposure map by outlining counties with high SVI (75th percentile or higher). Geographic unit will be county as this is the unit for glyphosate application rate. I will create separate maps for different crop types (i.e. soybeans, wheat, corn, cotton, oats, beans). If urban glyphosate usage data is available, I plan to use graduated symbols to highlight cities where its use is highest.
TO-DO
– does “No data” mean 0 or that they didn’t look at those states? – Look at the AHS air conditioning data: https://www.census.gov/programs-surveys/ahs/data/interactive/ahstablecreator.html?s_areas=42660&s_year=2021&s_tablename=TABLE3&s_bygroup1=1&s_bygroup2=1&s_filtergroup1=1&s_filtergroup2=1 – Figure out how to outline each county – Get the actual best heat measure that I want to use – Most prevalent type of gly resistant crop for county – Urban gly usage data?
Functions
plot_gly <- function(year, est_type) {
if (year==2013) {
data = merged_2013
}
if (year==2014) {
data = merged_2014
}
if (year==2015) {
data = merged_2015
}
if (year==2016) {
data = merged_2016
}
if (year==2017) {
data = merged_2017
}
if (est_type=="low") {
map <- ggplot(data) +
geom_sf(aes(fill = gly_per_land_low), linetype=0) +
ggtitle("Glyphosate per Land Area (Low Est.)", subtitle=year) +
theme_void(base_size = 14) +
theme(
plot.title = element_text(hjust = 0.5),
axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.ticks = element_blank(),
rect = element_blank()
) +
scale_fill_gradientn("Glyphosate usage \n(kg/km^2)",
colors = met.brewer("Tam"),
na.value = "grey50")
}
if (est_type=="high") {
map <- ggplot(data) +
geom_sf(aes(fill = gly_per_land_high), linetype=0) +
ggtitle("Glyphosate per Land Area (High Est.)", subtitle=year) +
theme_void(base_size = 14) +
theme(
plot.title = element_text(hjust = 0.5),
axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.ticks = element_blank(),
rect = element_blank()
) +
scale_fill_gradientn("Glyphosate usage \n(kg/km^2)",
colors = met.brewer("Tam"),
na.value = "grey50")
}
map
}
compounddata_formap <- function(compound, year) {
usgs_data <- usgs_ests %>% filter(COMPOUND==compound) %>%
filter(YEAR==year) %>% mutate(STATE_FIPS_CODE =
str_pad(as.character(STATE_FIPS_CODE),
2, pad = "0"),
COUNTY_FIPS_CODE =
str_pad(as.character(COUNTY_FIPS_CODE),
3, pad = "0"),
fips = as.factor(paste0(STATE_FIPS_CODE,
COUNTY_FIPS_CODE)))
if (year==2013) {
us_county_data <- us_county_2013
prism_county_data <- prism_county_2013_85s
svi_data <- svi_2014
}
if (year==2014) {
prism_county_data <- prism_county_2014_85s
svi_data <- svi_2014
}
if (year==2015) {
prism_county_data <- prism_county_2015_85s
svi_data <- svi_2016
}
if (year==2016) {
prism_county_data <- prism_county_2016_85s
svi_data <- svi_2016
}
if (year==2017) {
prism_county_data <- prism_county_2017_85s
svi_data <- svi_2018
}
merged_data <- left_join(us_county_2013, prism_county_data) %>%
left_join(usgs_data) %>% filter(STATEFP != "02" & STATEFP != "15" &
STATEFP != "72") %>%
mutate(epest_per_land_low = EPEST_LOW_KG/(ALAND/1e+6),
epest_per_land_high = EPEST_HIGH_KG/(ALAND/1e+6),
epest_per_totalarea_low = EPEST_LOW_KG/((ALAND+AWATER)/1e+6),
epest_per_totalarea_high = EPEST_HIGH_KG/((ALAND+AWATER)/1e+6)) %>%
left_join(svi_data)
merged_data
}
hot_compounddata_formap <- function(compound, year, temp_c) {
compounddata_formap(compound, year) %>% filter(tmean85_spray>=temp_c)
}
all_leaflet_low <- function(compound, year) {
#get data
data <- compounddata_formap(compound, year)
map <- data %>%
st_transform(crs = "+init=epsg:4326") %>%
leaflet(width = "100%") %>%
addProviderTiles(provider = "CartoDB.Voyager") %>%
addPolygons(
stroke = FALSE,
smoothFactor = 0,
fillOpacity = 0.5,
fillColor= ~
colorNumeric("YlOrRd", epest_per_land_low)(epest_per_land_low)) %>%
addLegend(pal = colorNumeric("YlOrRd", data$epest_per_land_low), # Define color palette
values = data$epest_per_land_low, # Values to map
title = "Pesticide Usage (kg/km^2)", # Legend title
position = "bottomright", # Position of the legend
na.label = "No data"
)
map
}
all_leaflet_high <- function(compound, year) {
#get data
data <- compounddata_formap(compound, year)
map <- data %>%
st_transform(crs = "+init=epsg:4326") %>%
leaflet(width = "100%") %>%
addProviderTiles(provider = "CartoDB.Voyager") %>%
addPolygons(
stroke = FALSE,
smoothFactor = 0,
fillOpacity = 0.5,
fillColor= ~
colorNumeric("YlOrRd", epest_per_land_high)(epest_per_land_high)) %>%
addLegend(pal = colorNumeric("YlOrRd", data$epest_per_land_high), # Define color palette
values = data$epest_per_land_high, # Values to map
title = "Pesticide Usage (kg/km^2)", # Legend title
position = "bottomright", # Position of the legend
na.label = "No data"
)
map
}
hot_leaflet_low <- function(compound, year, temp_c) {
#get data
data <- hot_compounddata_formap(compound, year, temp_c)
map <- data %>%
st_transform(crs = "+init=epsg:4326") %>%
leaflet(width = "100%") %>%
addProviderTiles(provider = "CartoDB.Voyager") %>%
addPolygons(
stroke = FALSE,
smoothFactor = 0,
fillOpacity = 0.5,
fillColor= ~
colorNumeric("YlOrRd", epest_per_land_low)(epest_per_land_low)) %>%
addLegend(pal = colorNumeric("YlOrRd", data$epest_per_land_low), # Define color palette
values = data$epest_per_land_low, # Values to map
title = "Pesticide Usage (kg/km^2)", # Legend title
position = "bottomright", # Position of the legend
na.label = "No data"
)
map
}
hot_leaflet_high <- function(compound, year, temp_c) {
#get data
data <- hot_compounddata_formap(compound, year, temp_c)
map <- data %>%
st_transform(crs = "+init=epsg:4326") %>%
leaflet(width = "100%") %>%
addProviderTiles(provider = "CartoDB.Voyager") %>%
addPolygons(
stroke = FALSE,
smoothFactor = 0,
fillOpacity = 0.5,
fillColor= ~
colorNumeric("YlOrRd", epest_per_land_high)(epest_per_land_high)) %>%
addLegend(pal = colorNumeric("YlOrRd", data$epest_per_land_high), # Define color palette
values = data$epest_per_land_high, # Values to map
title = "Pesticide Usage (kg/km^2)", # Legend title
position = "bottomright", # Position of the legend
na.label = "No data"
)
map
}Data Used
I am using USGS pesticide estimate data, Census Bureau data on U.S. county boundaries, and historical heat records in the U.S. from PRISM.
USGS Estimated Annual Agricultural Pesticide Use for Counties of the Conterminous United States, 2013-2017 (ver. 2.0, May 2020) data were obtained from here: https://www.sciencebase.gov/catalog/item/5e95c12282ce172707f2524e. Cartographic Boundary Files from the U.S Census Bureau were obtained from here for the year 2013: https://www.census.gov/geographies/mapping-files/time-series/geo/carto-boundary-file.2013.html#list-tab-1556094155. PRISM daily mean temperature data for 2013 was obtained from Robbie Parks’ Github (https://github.com/rmp15/PRISM-grids-into-FIPS-ZIP-censustract-USA/tree/main/output/fips/tmean). Citation for this work is Tuholske, C., Lynch, V.D., Spriggs, R. et al. Hazardous heat exposure among incarcerated people in the United States. Nat Sustain 7, 394–398 (2024). https://doi.org/10.1038/s41893-024-01293-y. Data dictionaries can be found at the above links.
I obtained CDC SVI data from: https://www.atsdr.cdc.gov/placeandhealth/svi/data_documentation_download.html. I obtained data for 2014, 2016, and 2018, and will use 2014 for years 2013-14, 2016 for 2015-16, and 2018 for 2017-2018. Data documentation is available at the same link (https://www.atsdr.cdc.gov/placeandhealth/svi/data_documentation_download.html).
#usgs_1992to2012 #read all these files from the list of files, then rbind them
usgs_ests <- fread(here("mapping data/EPest_county_estimates_2013_2017_v2.txt"))
usgs_compounds <- usgs_ests %>% select(c(COMPOUND)) %>% unique()
usgs_gly <- usgs_ests %>% filter(COMPOUND=="GLYPHOSATE")
usgs_gly_2013 <- usgs_gly %>% filter(YEAR==2013)
usgs_gly_2014 <- usgs_gly %>% filter(YEAR==2014)
usgs_gly_2015 <- usgs_gly %>% filter(YEAR==2015)
usgs_gly_2016 <- usgs_gly %>% filter(YEAR==2016)
usgs_gly_2017 <- usgs_gly %>% filter(YEAR==2017)
glimpse(usgs_gly)## Rows: 15,314
## Columns: 6
## $ COMPOUND <chr> "GLYPHOSATE", "GLYPHOSATE", "GLYPHOSATE", "GLYPHOSATE…
## $ YEAR <int> 2013, 2013, 2013, 2013, 2013, 2013, 2013, 2013, 2013,…
## $ STATE_FIPS_CODE <int> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,…
## $ COUNTY_FIPS_CODE <int> 1, 3, 5, 7, 9, 11, 13, 15, 17, 19, 21, 23, 25, 27, 29…
## $ EPEST_LOW_KG <dbl> 13336.7, 48158.8, 9683.4, 379.4, 5589.1, 7321.7, 2719…
## $ EPEST_HIGH_KG <dbl> 13538.5, 49014.8, 9683.4, 380.9, 5675.6, 7442.3, 2825…dim(usgs_gly)## [1] 15314 6us_county_2013 <- st_read(here("mapping data/cb_2013_us_county_20m/cb_2013_us_county_20m.shp"))## Reading layer `cb_2013_us_county_20m' from data source
## `/home/NETID/cm763/pesticide-heat-mapping/mapping data/cb_2013_us_county_20m/cb_2013_us_county_20m.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 3221 features and 9 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -179.1473 ymin: 17.88481 xmax: 179.7785 ymax: 71.35256
## Geodetic CRS: NAD83glimpse(us_county_2013)## Rows: 3,221
## Columns: 10
## $ STATEFP <chr> "06", "13", "18", "21", "21", "26", "01", "48", "01", "36", "…
## $ COUNTYFP <chr> "049", "257", "127", "077", "053", "149", "085", "485", "037"…
## $ COUNTYNS <chr> "00277289", "00350028", "00450382", "00516885", "00516873", "…
## $ AFFGEOID <chr> "0500000US06049", "0500000US13257", "0500000US18127", "050000…
## $ GEOID <chr> "06049", "13257", "18127", "21077", "21053", "26149", "01085"…
## $ NAME <chr> "Modoc", "Stephens", "Porter", "Gallatin", "Clinton", "St. Jo…
## $ LSAD <chr> "06", "06", "06", "06", "06", "06", "06", "06", "06", "06", "…
## $ ALAND <dbl> 10140841147, 463948108, 1083011641, 254698298, 510864239, 129…
## $ AWATER <dbl> 745929078, 13121937, 268393387, 16519758, 21164150, 52934828,…
## $ geometry <MULTIPOLYGON [°]> MULTIPOLYGON (((-121.4475 4..., MULTIPOLYGON (((…dim(us_county_2013)## [1] 3221 10prism_county_2013 <- readRDS(here("mapping data/weighted_area_raster_fips_tmean_daily_2013.rds"))
glimpse(prism_county_2013)## Rows: 1,134,420
## Columns: 6
## $ fips <fct> 01001, 01003, 01005, 01007, 01009, 01011, 01013, 01015, 01017, 0…
## $ tmean <dbl> 6.188720, 9.711895, 7.783865, 6.224829, 5.372772, 6.969209, 7.23…
## $ date <chr> "01/01/2013", "01/01/2013", "01/01/2013", "01/01/2013", "01/01/2…
## $ day <chr> "01", "01", "01", "01", "01", "01", "01", "01", "01", "01", "01"…
## $ month <chr> "01", "01", "01", "01", "01", "01", "01", "01", "01", "01", "01"…
## $ year <dbl> 2013, 2013, 2013, 2013, 2013, 2013, 2013, 2013, 2013, 2013, 2013…dim(prism_county_2013)## [1] 1134420 6prism_county_2014 <- readRDS(here("mapping data/weighted_area_raster_fips_tmean_daily_2014.rds"))
prism_county_2015 <- readRDS(here("mapping data/weighted_area_raster_fips_tmean_daily_2015.rds"))
prism_county_2016 <- readRDS(here("mapping data/weighted_area_raster_fips_tmean_daily_2016.rds"))
prism_county_2017 <- readRDS(here("mapping data/weighted_area_raster_fips_tmean_daily_2017.rds"))
svi_2014 <- fread(here("mapping data/SVI_2014_US_county.csv"))
glimpse(svi_2014)## Rows: 3,142
## Columns: 127
## $ AFFGEOID <chr> "0500000US01001", "0500000US01003", "0500000US01005…
## $ ST <int> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, …
## $ STATE <chr> "ALABAMA", "ALABAMA", "ALABAMA", "ALABAMA", "ALABAM…
## $ ST_ABBR <chr> "AL", "AL", "AL", "AL", "AL", "AL", "AL", "AL", "AL…
## $ COUNTY <chr> "Autauga", "Baldwin", "Barbour", "Bibb", "Blount", …
## $ FIPS <int> 1001, 1003, 1005, 1007, 1009, 1011, 1013, 1015, 101…
## $ LOCATION <chr> "Autauga County, Alabama", "Baldwin County, Alabama…
## $ AREA_SQMI <dbl> 594.4366, 1589.8073, 884.8767, 622.5824, 644.8065, …
## $ E_TOTPOP <dbl> 55136, 191205, 27119, 22653, 57645, 10693, 20523, 1…
## $ M_TOTPOP <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ E_HU <dbl> 22431, 105563, 11833, 8985, 23868, 4469, 9934, 5330…
## $ M_HU <dbl> 67, 168, 120, 66, 77, 100, 67, 209, 50, 89, 68, 33,…
## $ E_HH <dbl> 20304, 73058, 9145, 7078, 20934, 3746, 8253, 45348,…
## $ M_HH <dbl> 458, 1241, 311, 390, 399, 216, 252, 705, 393, 451, …
## $ E_POV <dbl> 7006, 25988, 5832, 3596, 9866, 2085, 5239, 24794, 8…
## $ M_POV <dbl> 935, 2457, 639, 770, 947, 477, 517, 1491, 823, 876,…
## $ E_UNEMP <dbl> 2252, 7856, 1527, 975, 2291, 809, 1075, 7257, 1916,…
## $ M_UNEMP <dbl> 352, 918, 282, 310, 358, 187, 177, 594, 304, 221, 3…
## $ E_PCI <dbl> 24644, 26851, 17350, 18110, 20501, 17706, 18115, 21…
## $ M_PCI <dbl> 780, 712, 821, 1477, 719, 1557, 832, 573, 1482, 138…
## $ E_NOHSDP <dbl> 5012, 14615, 4790, 3466, 8567, 2511, 3288, 15674, 5…
## $ M_NOHSDP <dbl> 561, 1051, 341, 500, 697, 319, 278, 787, 388, 373, …
## $ E_AGE65 <dbl> 7321, 33782, 4180, 3209, 9172, 1529, 3592, 17820, 6…
## $ M_AGE65 <dbl> 84, 151, 81, 84, 89, 0, 72, 95, 38, 94, 135, 0, 0, …
## $ E_AGE17 <dbl> 14214, 43186, 5862, 4903, 13877, 2283, 4835, 26337,…
## $ M_AGE17 <dbl> 0, 0, 27, 0, 33, 0, 50, 0, 26, 78, 50, 0, 0, 46, 10…
## $ E_DISABL <dbl> 8700, 26603, 4856, 3343, 9729, 2049, 4301, 21597, 6…
## $ M_DISABL <dbl> 674, 1468, 379, 429, 748, 309, 284, 861, 507, 552, …
## $ E_SNGPNT <dbl> 1649, 5027, 1098, 637, 1650, 535, 1104, 4406, 1725,…
## $ M_SNGPNT <dbl> 300.2, 603.7, 160.7, 211.3, 299.2, 179.6, 173.7, 50…
## $ E_MINRTY <dbl> 13103, 32158, 14614, 5717, 6825, 8321, 9517, 31315,…
## $ M_MINRTY <dbl> 109, 162, 166, 21, 52, 18, 22, 68, 24, 21, 91, 18, …
## $ E_LIMENG <dbl> 249, 2571, 549, 112, 954, 257, 82, 915, 67, 25, 113…
## $ M_LIMENG <dbl> 138.5, 658.6, 261.6, 110.9, 222.2, 250.7, 100.2, 24…
## $ E_MUNIT <dbl> 842, 18988, 129, 163, 215, 65, 136, 2025, 755, 195,…
## $ M_MUNIT <dbl> 344.4, 1015.0, 68.6, 97.6, 95.2, 50.6, 67.7, 296.9,…
## $ E_MOBILE <dbl> 4401, 12200, 3175, 2379, 5497, 1438, 2397, 7467, 30…
## $ M_MOBILE <dbl> 367, 821, 270, 314, 482, 228, 203, 511, 320, 398, 4…
## $ E_CROWD <dbl> 530, 998, 178, 17, 344, 122, 95, 575, 506, 266, 567…
## $ M_CROWD <dbl> 231.2, 265.6, 83.6, 25.8, 142.4, 118.7, 49.2, 165.0…
## $ E_NOVEH <dbl> 1081, 2242, 802, 299, 823, 703, 701, 2562, 1238, 34…
## $ M_NOVEH <dbl> 218, 308, 152, 103, 180, 188, 127, 307, 220, 98, 18…
## $ E_GROUPQ <dbl> 442, 2611, 2869, 1576, 572, 576, 270, 2768, 447, 39…
## $ M_GROUPQ <dbl> 133, 405, 290, 177, 147, 248, 146, 500, 105, 187, 1…
## $ EP_POV <dbl> 12.8, 13.8, 24.1, 17.0, 17.3, 20.5, 25.9, 21.7, 23.…
## $ MP_POV <dbl> 1.7, 1.3, 2.6, 3.6, 1.7, 4.7, 2.6, 1.3, 2.4, 3.4, 2…
## $ EP_UNEMP <dbl> 8.5, 8.6, 14.2, 10.9, 9.3, 17.4, 12.2, 13.5, 12.7, …
## $ MP_UNEMP <dbl> 1.3, 1.0, 2.5, 3.4, 1.4, 3.6, 1.9, 1.1, 2.0, 2.0, 1…
## $ EP_PCI <dbl> 24644, 26851, 17350, 18110, 20501, 17706, 18115, 21…
## $ MP_PCI <dbl> 780, 712, 821, 1477, 719, 1557, 832, 573, 1482, 138…
## $ EP_NOHSDP <dbl> 13.8, 11.0, 25.4, 22.1, 21.9, 34.5, 23.4, 19.9, 22.…
## $ MP_NOHSDP <dbl> 1.5, 0.8, 1.8, 3.2, 1.8, 4.4, 2.0, 1.0, 1.6, 2.0, 1…
## $ EP_AGE65 <dbl> 13.3, 17.7, 15.4, 14.2, 15.9, 14.3, 17.5, 15.2, 17.…
## $ MP_AGE65 <dbl> 0.2, 0.1, 0.3, 0.4, 0.2, 0.0, 0.3, 0.1, 0.1, 0.4, 0…
## $ EP_AGE17 <dbl> 25.8, 22.6, 21.6, 21.6, 24.1, 21.4, 23.6, 22.5, 21.…
## $ MP_AGE17 <dbl> 0.0, 0.0, 0.1, 0.0, 0.1, 0.0, 0.2, 0.0, 0.1, 0.3, 0…
## $ EP_DISABL <dbl> 16.0, 14.1, 20.0, 15.8, 17.0, 20.2, 21.2, 18.7, 19.…
## $ MP_DISABL <dbl> 1.2, 0.8, 1.5, 2.0, 1.3, 3.0, 1.4, 0.7, 1.5, 2.2, 1…
## $ EP_SNGPNT <dbl> 8.1, 6.9, 12.0, 9.0, 7.9, 14.3, 13.4, 9.7, 12.4, 8.…
## $ MP_SNGPNT <dbl> 1.5, 0.8, 1.7, 2.9, 1.4, 4.7, 2.1, 1.1, 1.6, 2.0, 1…
## $ EP_MINRTY <dbl> 23.8, 16.8, 53.9, 25.2, 11.8, 77.8, 46.4, 26.7, 42.…
## $ MP_MINRTY <dbl> 0.2, 0.1, 0.6, 0.1, 0.1, 0.2, 0.1, 0.1, 0.1, 0.1, 0…
## $ EP_LIMENG <dbl> 0.5, 1.4, 2.2, 0.5, 1.8, 2.6, 0.4, 0.8, 0.2, 0.1, 2…
## $ MP_LIMENG <dbl> 0.3, 0.4, 1.0, 0.5, 0.4, 2.5, 0.5, 0.2, 0.3, 0.3, 0…
## $ EP_MUNIT <dbl> 3.8, 18.0, 1.1, 1.8, 0.9, 1.5, 1.4, 3.8, 4.5, 1.2, …
## $ MP_MUNIT <dbl> 1.5, 1.0, 0.6, 1.1, 0.4, 1.1, 0.7, 0.6, 1.0, 0.5, 0…
## $ EP_MOBILE <dbl> 19.6, 11.6, 26.8, 26.5, 23.0, 32.2, 24.1, 14.0, 18.…
## $ MP_MOBILE <dbl> 1.6, 0.8, 2.3, 3.5, 2.0, 5.0, 2.1, 1.0, 1.9, 2.5, 2…
## $ EP_CROWD <dbl> 2.6, 1.4, 1.9, 0.2, 1.6, 3.3, 1.2, 1.3, 3.6, 2.3, 3…
## $ MP_CROWD <dbl> 1.1, 0.4, 0.9, 0.4, 0.7, 3.2, 0.6, 0.4, 1.0, 1.2, 1…
## $ EP_NOVEH <dbl> 5.3, 3.1, 8.8, 4.2, 3.9, 18.8, 8.5, 5.6, 8.9, 3.0, …
## $ MP_NOVEH <dbl> 1.1, 0.4, 1.6, 1.5, 0.9, 5.0, 1.5, 0.7, 1.6, 0.8, 1…
## $ EP_GROUPQ <dbl> 0.8, 1.4, 10.6, 7.0, 1.0, 5.4, 1.3, 2.4, 1.3, 1.5, …
## $ MP_GROUPQ <dbl> 0.2, 0.2, 1.1, 0.8, 0.3, 2.3, 0.7, 0.4, 0.3, 0.7, 0…
## $ EPL_POV <dbl> 0.2833, 0.3531, 0.8736, 0.5565, 0.5782, 0.7561, 0.9…
## $ EPL_UNEMP <dbl> 0.5174, 0.5288, 0.9319, 0.7650, 0.6074, 0.9787, 0.8…
## $ EPL_PCI <dbl> 0.3913, 0.2467, 0.9274, 0.8949, 0.7348, 0.9115, 0.8…
## $ EPL_NOHSDP <dbl> 0.5126, 0.3225, 0.9236, 0.8360, 0.8287, 0.9901, 0.8…
## $ SPL_THEME1 <dbl> 1.7046, 1.4511, 3.6565, 3.0525, 2.7491, 3.6364, 3.5…
## $ RPL_THEME1 <dbl> 0.4145, 0.3372, 0.9628, 0.8201, 0.7348, 0.9580, 0.9…
## $ EPL_AGE65 <dbl> 0.1980, 0.6281, 0.3964, 0.2677, 0.4429, 0.2767, 0.6…
## $ EPL_AGE17 <dbl> 0.8532, 0.4597, 0.3327, 0.3368, 0.6714, 0.3044, 0.6…
## $ EPL_DISABL <dbl> 0.5762, 0.3983, 0.8475, 0.5571, 0.6587, 0.8580, 0.8…
## $ EPL_SNGPNT <dbl> 0.4330, 0.2372, 0.8937, 0.5804, 0.3903, 0.9659, 0.9…
## $ SPL_THEME2 <dbl> 2.0605, 1.7233, 2.4702, 1.7421, 2.1633, 2.4050, 3.0…
## $ RPL_THEME2 <dbl> 0.5387, 0.2811, 0.8195, 0.2948, 0.6199, 0.7826, 0.9…
## $ EPL_MINRTY <dbl> 0.6332, 0.5336, 0.9086, 0.6517, 0.4247, 0.9819, 0.8…
## $ EPL_LIMENG <dbl> 0.3722, 0.6778, 0.7720, 0.3961, 0.7297, 0.8042, 0.3…
## $ SPL_THEME3 <dbl> 1.0054, 1.2114, 1.6807, 1.0478, 1.1544, 1.7861, 1.2…
## $ RPL_THEME3 <dbl> 0.4986, 0.6288, 0.8898, 0.5237, 0.5864, 0.9322, 0.6…
## $ EPL_MUNIT <dbl> 0.6135, 0.9733, 0.2006, 0.3499, 0.1649, 0.2767, 0.2…
## $ EPL_MOBILE <dbl> 0.7701, 0.5183, 0.9064, 0.9019, 0.8488, 0.9577, 0.8…
## $ EPL_CROWD <dbl> 0.7011, 0.3098, 0.5317, 0.0216, 0.4167, 0.8103, 0.2…
## $ EPL_NOVEH <dbl> 0.3820, 0.0844, 0.8341, 0.2073, 0.1713, 0.9898, 0.8…
## $ EPL_GROUPQ <dbl> 0.0863, 0.2865, 0.9414, 0.8781, 0.1407, 0.8259, 0.2…
## $ SPL_THEME4 <dbl> 2.5530, 2.1722, 3.4142, 2.3588, 1.7424, 3.8602, 2.4…
## $ RPL_THEME4 <dbl> 0.5100, 0.3238, 0.9089, 0.4152, 0.1563, 0.9866, 0.4…
## $ SPL_THEMES <dbl> 7.3235, 6.5581, 11.2216, 8.2012, 7.8093, 11.6877, 1…
## $ RPL_THEMES <dbl> 0.4696, 0.3432, 0.9742, 0.6278, 0.5581, 0.9914, 0.9…
## $ F_POV <dbl> 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, …
## $ F_UNEMP <dbl> 0, 0, 1, 0, 0, 1, 0, 1, 0, 0, 0, 1, 1, 0, 0, 0, 0, …
## $ F_PCI <dbl> 0, 0, 1, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ F_NOHSDP <dbl> 0, 0, 1, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 1, 1, 0, 0, …
## $ F_THEME1 <dbl> 0, 0, 3, 0, 0, 3, 1, 1, 0, 0, 0, 2, 2, 1, 1, 0, 0, …
## $ F_AGE65 <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ F_AGE17 <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ F_DISABL <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 1, 0, 0, 0, 0, 0, …
## $ F_SNGPNT <dbl> 0, 0, 0, 0, 0, 1, 1, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ F_THEME2 <dbl> 0, 0, 0, 0, 0, 1, 1, 0, 1, 1, 0, 1, 0, 0, 0, 0, 0, …
## $ F_MINRTY <dbl> 0, 0, 1, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ F_LIMENG <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ F_THEME3 <dbl> 0, 0, 1, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ F_MUNIT <dbl> 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ F_MOBILE <dbl> 0, 0, 1, 1, 0, 1, 0, 0, 0, 1, 1, 1, 0, 1, 1, 0, 0, …
## $ F_CROWD <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ F_NOVEH <dbl> 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, …
## $ F_GROUPQ <dbl> 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ F_THEME4 <dbl> 0, 1, 2, 1, 0, 2, 0, 0, 0, 1, 1, 1, 1, 1, 1, 0, 0, …
## $ F_TOTAL <dbl> 0, 1, 6, 1, 0, 7, 2, 1, 1, 2, 1, 4, 3, 2, 2, 0, 0, …
## $ E_UNINSUR <dbl> 5814, 27758, 3760, 2725, 6922, 1791, 2603, 15258, 5…
## $ M_UNINSUR <dbl> 674, 2097, 412, 554, 660, 313, 327, 1031, 476, 504,…
## $ EP_UNINSUR <dbl> 10.7, 14.7, 15.5, 12.9, 12.1, 17.6, 12.8, 13.2, 15.…
## $ MP_UNINSUR <dbl> 1.2, 1.1, 1.7, 2.6, 1.2, 3.0, 1.6, 0.9, 1.4, 2.0, 1…
## $ E_DAYPOP <dbl> 43534, 177010, 29769, 19274, 41857, 10639, 19896, 1…
## $ Shape <chr> "(-86.64274304015402, 32.53492165233667)", "(-87.72…
## $ `Shape.STArea()` <dbl> 0.1502601, 0.4099377, 0.2232557, 0.1565251, 0.16440…
## $ `Shape.STLength()` <dbl> 2.052620, 4.277985, 2.566611, 1.886957, 2.392305, 2…dim(svi_2014)## [1] 3142 127svi_2016 <- fread(here("mapping data/SVI_2016_US_county.csv"))
glimpse(svi_2016)## Rows: 3,142
## Columns: 123
## $ ST <int> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,…
## $ STATE <chr> "ALABAMA", "ALABAMA", "ALABAMA", "ALABAMA", "ALABAMA", "ALA…
## $ ST_ABBR <chr> "AL", "AL", "AL", "AL", "AL", "AL", "AL", "AL", "AL", "AL",…
## $ COUNTY <chr> "Autauga", "Blount", "Chambers", "Coffee", "Colbert", "Covi…
## $ FIPS <int> 1001, 1009, 1017, 1031, 1033, 1039, 1043, 1045, 1051, 1055,…
## $ LOCATION <chr> "Autauga County, Alabama", "Blount County, Alabama", "Chamb…
## $ AREA_SQMI <dbl> 594.4461, 644.8065, 596.5311, 678.9857, 592.6197, 1030.4558…
## $ E_TOTPOP <dbl> 55049, 57704, 34018, 50991, 54377, 37729, 81316, 49607, 809…
## $ M_TOTPOP <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ E_HU <dbl> 22714, 23850, 16905, 22862, 26156, 18783, 37150, 22839, 333…
## $ M_HU <dbl> 75, 59, 53, 106, 88, 66, 92, 95, 217, 126, 63, 138, 63, 146…
## $ E_HH <dbl> 20800, 20619, 13851, 19375, 22105, 15179, 31081, 18794, 288…
## $ M_HH <dbl> 391, 403, 359, 446, 428, 337, 517, 416, 541, 505, 246, 369,…
## $ E_POV <dbl> 6697, 9441, 6805, 8219, 8910, 7031, 13718, 9897, 10893, 187…
## $ M_POV <dbl> 1037, 963, 849, 901, 839, 644, 1340, 914, 1217, 1295, 455, …
## $ E_UNEMP <dbl> 1437, 1367, 1136, 1410, 1795, 1723, 2081, 2145, 2694, 3837,…
## $ M_UNEMP <dbl> 277, 284, 232, 235, 301, 314, 314, 302, 522, 462, 170, 224,…
## $ E_PCI <dbl> 26168, 21033, 21532, 25325, 23318, 21738, 21041, 22834, 247…
## $ M_PCI <dbl> 1221, 689, 1399, 939, 854, 826, 734, 783, 1013, 489, 1459, …
## $ E_NOHSDP <dbl> 4528, 7882, 4708, 5145, 6344, 5051, 10050, 4726, 7238, 1258…
## $ M_NOHSDP <dbl> 445, 645, 421, 418, 506, 345, 745, 396, 594, 619, 239, 408,…
## $ E_AGE65 <dbl> 7695, 9921, 6255, 8048, 10034, 7613, 14435, 7654, 11243, 17…
## $ M_AGE65 <dbl> 104, 123, 48, 89, 71, 69, 106, 79, 82, 75, 46, 0, 29, 78, 5…
## $ E_AGE17 <dbl> 13853, 13601, 7283, 12122, 11735, 8318, 18279, 11728, 18553…
## $ M_AGE17 <dbl> 34, 29, 52, 0, 0, 67, 110, 0, 0, 107, 47, 0, 66, 24, 60, 63…
## $ E_DISABL <dbl> 10009, 8538, 5960, 8942, 10561, 7858, 12621, 9865, 13201, 2…
## $ M_DISABL <dbl> 850, 663, 492, 597, 677, 439, 739, 494, 967, 937, 296, 462,…
## $ E_SNGPNT <dbl> 1516, 1614, 1354, 2018, 1879, 1240, 2168, 1649, 2688, 3631,…
## $ M_SNGPNT <dbl> 267.0, 301.0, 187.6, 266.6, 290.0, 171.4, 337.1, 251.6, 394…
## $ E_MINRTY <dbl> 13386, 7122, 14715, 14598, 11499, 6220, 6128, 14787, 21301,…
## $ M_MINRTY <dbl> 161, 147, 24, 7, 46, 6, 104, 77, 109, 108, 18, 55, 91, 58, …
## $ E_LIMENG <dbl> 432, 1018, 114, 716, 104, 211, 746, 395, 635, 1475, 115, 73…
## $ M_LIMENG <dbl> 163.3, 248.1, 100.7, 242.1, 122.2, 112.3, 207.7, 183.9, 247…
## $ E_MUNIT <dbl> 1034, 190, 597, 218, 641, 174, 901, 304, 894, 1515, 86, 199…
## $ M_MUNIT <dbl> 329.9, 85.8, 170.3, 81.8, 154.7, 75.9, 191.1, 113.8, 235.8,…
## $ E_MOBILE <dbl> 4095, 5467, 2695, 2863, 2478, 3505, 7996, 4020, 5587, 5487,…
## $ M_MOBILE <dbl> 379, 429, 317, 202, 333, 255, 575, 302, 472, 441, 239, 297,…
## $ E_CROWD <dbl> 254, 391, 555, 270, 150, 237, 547, 269, 332, 639, 59, 479, …
## $ M_CROWD <dbl> 104.5, 130.6, 150.8, 101.0, 76.9, 90.7, 136.0, 101.2, 129.8…
## $ E_NOVEH <dbl> 1024, 816, 1110, 1126, 1361, 797, 1264, 1198, 1291, 2513, 4…
## $ M_NOVEH <dbl> 242, 198, 183, 204, 224, 123, 218, 213, 308, 286, 117, 155,…
## $ E_GROUPQ <dbl> 490, 552, 482, 616, 432, 612, 1164, 1209, 4651, 1403, 268, …
## $ M_GROUPQ <dbl> 163, 131, 106, 134, 102, 151, 216, 219, 570, 249, 79, 170, …
## $ EP_POV <dbl> 12.3, 16.5, 20.3, 16.4, 16.5, 19.0, 17.1, 20.5, 14.3, 18.4,…
## $ MP_POV <dbl> 1.9, 1.7, 2.5, 1.8, 1.6, 1.7, 1.7, 1.9, 1.6, 1.3, 2.7, 0.7,…
## $ EP_UNEMP <dbl> 5.6, 6.0, 7.5, 6.2, 7.5, 10.7, 6.1, 10.4, 7.5, 8.5, 7.5, 8.…
## $ MP_UNEMP <dbl> 1.1, 1.2, 1.5, 1.0, 1.3, 1.9, 0.9, 1.4, 1.4, 1.0, 2.2, 0.5,…
## $ EP_PCI <dbl> 26168, 21033, 21532, 25325, 23318, 21738, 21041, 22834, 247…
## $ MP_PCI <dbl> 1221, 689, 1399, 939, 854, 826, 734, 783, 1013, 489, 1459, …
## $ EP_NOHSDP <dbl> 12.4, 20.0, 19.7, 14.9, 16.6, 19.1, 17.8, 14.2, 13.2, 17.5,…
## $ MP_NOHSDP <dbl> 1.2, 1.6, 1.8, 1.2, 1.3, 1.3, 1.3, 1.2, 1.1, 0.9, 1.9, 0.6,…
## $ EP_AGE65 <dbl> 14.0, 17.2, 18.4, 15.8, 18.5, 20.2, 17.8, 15.4, 13.9, 17.4,…
## $ MP_AGE65 <dbl> 0.2, 0.2, 0.1, 0.2, 0.1, 0.2, 0.1, 0.2, 0.1, 0.1, 0.3, 0.0,…
## $ EP_AGE17 <dbl> 25.2, 23.6, 21.4, 23.8, 21.6, 22.0, 22.5, 23.6, 22.9, 22.0,…
## $ MP_AGE17 <dbl> 0.1, 0.1, 0.2, 0.0, 0.0, 0.2, 0.1, 0.0, 0.0, 0.1, 0.3, 0.0,…
## $ EP_DISABL <dbl> 18.4, 14.9, 17.7, 18.2, 19.6, 21.1, 15.7, 21.0, 17.4, 19.8,…
## $ MP_DISABL <dbl> 1.6, 1.2, 1.5, 1.2, 1.3, 1.2, 0.9, 1.1, 1.3, 0.9, 1.7, 0.4,…
## $ EP_SNGPNT <dbl> 7.3, 7.8, 9.8, 10.4, 8.5, 8.2, 7.0, 8.8, 9.3, 9.2, 7.2, 11.…
## $ MP_SNGPNT <dbl> 1.3, 1.5, 1.3, 1.4, 1.3, 1.1, 1.1, 1.3, 1.4, 0.9, 1.9, 0.6,…
## $ EP_MINRTY <dbl> 24.3, 12.3, 43.3, 28.6, 21.1, 16.5, 7.5, 29.8, 26.3, 21.6, …
## $ MP_MINRTY <dbl> 0.3, 0.3, 0.1, 0.0, 0.1, 0.0, 0.1, 0.2, 0.1, 0.1, 0.1, 0.1,…
## $ EP_LIMENG <dbl> 0.8, 1.9, 0.4, 1.5, 0.2, 0.6, 1.0, 0.9, 0.8, 1.5, 0.7, 0.8,…
## $ MP_LIMENG <dbl> 0.3, 0.5, 0.3, 0.5, 0.2, 0.3, 0.3, 0.4, 0.3, 0.3, 0.5, 0.1,…
## $ EP_MUNIT <dbl> 4.6, 0.8, 3.5, 1.0, 2.5, 0.9, 2.4, 1.3, 2.7, 3.2, 1.0, 4.3,…
## $ MP_MUNIT <dbl> 1.5, 0.4, 1.0, 0.4, 0.6, 0.4, 0.5, 0.5, 0.7, 0.4, 0.5, 0.4,…
## $ EP_MOBILE <dbl> 18.0, 22.9, 15.9, 12.5, 9.5, 18.7, 21.5, 17.6, 16.8, 11.6, …
## $ MP_MOBILE <dbl> 1.7, 1.8, 1.9, 0.9, 1.3, 1.4, 1.5, 1.3, 1.4, 0.9, 2.7, 0.6,…
## $ EP_CROWD <dbl> 1.2, 1.9, 4.0, 1.4, 0.7, 1.6, 1.8, 1.4, 1.1, 1.6, 0.9, 1.2,…
## $ MP_CROWD <dbl> 0.5, 0.6, 1.1, 0.5, 0.3, 0.6, 0.4, 0.5, 0.4, 0.5, 0.5, 0.2,…
## $ EP_NOVEH <dbl> 4.9, 4.0, 8.0, 5.8, 6.2, 5.3, 4.1, 6.4, 4.5, 6.4, 6.7, 6.3,…
## $ MP_NOVEH <dbl> 1.1, 1.0, 1.3, 1.0, 1.0, 0.8, 0.7, 1.1, 1.1, 0.7, 1.7, 0.4,…
## $ EP_GROUPQ <dbl> 0.9, 1.0, 1.4, 1.2, 0.8, 1.6, 1.4, 2.4, 5.7, 1.4, 1.6, 0.9,…
## $ MP_GROUPQ <dbl> 0.3, 0.2, 0.3, 0.3, 0.2, 0.4, 0.3, 0.4, 0.7, 0.2, 0.5, 0.2,…
## $ EPL_POV <dbl> 0.2824, 0.5536, 0.7644, 0.5466, 0.5536, 0.7042, 0.5925, 0.7…
## $ EPL_UNEMP <dbl> 0.3298, 0.3792, 0.5947, 0.4094, 0.5947, 0.8854, 0.3964, 0.8…
## $ EPL_PCI <dbl> 0.3607, 0.7504, 0.7122, 0.4187, 0.5766, 0.6969, 0.7494, 0.6…
## $ EPL_NOHSDP <dbl> 0.4744, 0.8058, 0.7934, 0.6055, 0.6829, 0.7720, 0.7329, 0.5…
## $ SPL_THEME1 <dbl> 1.4473, 2.4890, 2.8647, 1.9803, 2.4078, 3.0586, 2.4712, 2.8…
## $ RPL_THEME1 <dbl> 0.3298, 0.6609, 0.7689, 0.5072, 0.6345, 0.8214, 0.6558, 0.7…
## $ EPL_AGE65 <dbl> 0.1964, 0.4909, 0.6151, 0.3464, 0.6243, 0.7676, 0.5524, 0.3…
## $ EPL_AGE17 <dbl> 0.8313, 0.6466, 0.3467, 0.6756, 0.3664, 0.4263, 0.4938, 0.6…
## $ EPL_DISABL <dbl> 0.7380, 0.4527, 0.6890, 0.7233, 0.8032, 0.8768, 0.5288, 0.8…
## $ EPL_SNGPNT <dbl> 0.3200, 0.4018, 0.7138, 0.7845, 0.5244, 0.4629, 0.2744, 0.5…
## $ SPL_THEME2 <dbl> 2.0856, 1.9920, 2.3645, 2.5298, 2.3184, 2.5336, 1.8494, 2.4…
## $ RPL_THEME2 <dbl> 0.5568, 0.4664, 0.7615, 0.8507, 0.7326, 0.8516, 0.3613, 0.7…
## $ EPL_MINRTY <dbl> 0.6339, 0.4238, 0.8402, 0.6877, 0.5880, 0.5205, 0.2786, 0.6…
## $ EPL_LIMENG <dbl> 0.5355, 0.7482, 0.2996, 0.6931, 0.1821, 0.4323, 0.5813, 0.5…
## $ SPL_THEME3 <dbl> 1.1694, 1.1719, 1.1398, 1.3808, 0.7701, 0.9529, 0.8599, 1.2…
## $ RPL_THEME3 <dbl> 0.5976, 0.5992, 0.5798, 0.7313, 0.3642, 0.4686, 0.4174, 0.6…
## $ EPL_MUNIT <dbl> 0.6791, 0.1344, 0.5845, 0.1668, 0.4511, 0.1576, 0.4495, 0.2…
## $ EPL_MOBILE <dbl> 0.7268, 0.8465, 0.6734, 0.5613, 0.4368, 0.7488, 0.8138, 0.7…
## $ EPL_CROWD <dbl> 0.2477, 0.5056, 0.8812, 0.3091, 0.0780, 0.3798, 0.4550, 0.3…
## $ EPL_NOVEH <dbl> 0.3298, 0.1907, 0.7740, 0.4877, 0.5495, 0.3986, 0.2038, 0.5…
## $ EPL_GROUPQ <dbl> 0.1251, 0.1515, 0.3251, 0.2340, 0.0981, 0.4024, 0.3308, 0.5…
## $ SPL_THEME4 <dbl> 2.1086, 1.8287, 3.2381, 1.7590, 1.6135, 2.0872, 2.2528, 2.4…
## $ RPL_THEME4 <dbl> 0.2881, 0.1827, 0.8596, 0.1624, 0.1175, 0.2802, 0.3531, 0.4…
## $ SPL_THEMES <dbl> 6.8109, 7.4817, 9.6071, 7.6498, 7.1098, 8.6323, 7.4333, 8.9…
## $ RPL_THEMES <dbl> 0.3773, 0.4986, 0.8408, 0.5266, 0.4333, 0.6950, 0.4909, 0.7…
## $ F_POV <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_UNEMP <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_PCI <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_NOHSDP <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_THEME1 <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_AGE65 <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_AGE17 <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_DISABL <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_SNGPNT <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_THEME2 <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_MINRTY <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_LIMENG <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_THEME3 <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_MUNIT <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_MOBILE <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_CROWD <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_NOVEH <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_GROUPQ <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_THEME4 <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_TOTAL <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ E_UNINSUR <dbl> 4852, 6388, 3979, 5253, 4932, 4617, 11704, 6130, 6665, 1236…
## $ M_UNINSUR <dbl> 649, 740, 544, 464, 458, 448, 972, 562, 803, 857, 271, 597,…
## $ EP_UNINSUR <dbl> 8.9, 11.2, 11.8, 10.7, 9.1, 12.4, 14.5, 13.1, 8.8, 12.1, 11…
## $ MP_UNINSUR <dbl> 1.2, 1.3, 1.6, 0.9, 0.8, 1.2, 1.2, 1.2, 1.0, 0.8, 1.6, 0.6,…
## $ E_DAYPOP <dbl> 40854, 42597, 27940, 47236, 56227, 37246, 79900, 49165, 669…dim(svi_2016)## [1] 3142 123svi_2018 <- fread(here("mapping data/SVI_2018_US_county.csv"))
glimpse(svi_2018)## Rows: 3,142
## Columns: 123
## $ ST <int> 35, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1…
## $ STATE <chr> "NEW MEXICO", "ALABAMA", "ALABAMA", "ALABAMA", "ALABAMA", "…
## $ ST_ABBR <chr> "NM", "AL", "AL", "AL", "AL", "AL", "AL", "AL", "AL", "AL",…
## $ COUNTY <chr> "Rio Arriba", "Autauga", "Blount", "Butler", "Calhoun", "Ch…
## $ FIPS <int> 35039, 1001, 1009, 1013, 1015, 1017, 1031, 1033, 1039, 1043…
## $ LOCATION <chr> "Rio Arriba County, New Mexico", "Autauga County, Alabama",…
## $ AREA_SQMI <dbl> 5860.8692, 594.4435, 644.8305, 776.8382, 605.8673, 596.5606…
## $ E_TOTPOP <int> 39307, 55200, 57645, 20025, 115098, 33826, 51288, 54495, 37…
## $ M_TOTPOP <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ E_HU <int> 20044, 23315, 24222, 10026, 53682, 16981, 23088, 26466, 189…
## $ M_HU <int> 71, 71, 55, 51, 184, 72, 57, 58, 85, 92, 110, 192, 75, 44, …
## $ E_HH <int> 12398, 21115, 20600, 6708, 45033, 13516, 19789, 21799, 1500…
## $ M_HH <int> 439, 383, 396, 274, 683, 372, 339, 445, 381, 592, 402, 638,…
## $ E_POV <int> -999, 8422, 8220, 4640, 20819, 5531, 7604, 8785, 6702, 1283…
## $ M_POV <int> -999, 1137, 992, 521, 1317, 682, 731, 848, 683, 1294, 859, …
## $ E_UNEMP <int> -999, 1065, 909, 567, 4628, 773, 1371, 1397, 1428, 1903, 17…
## $ M_UNEMP <int> -999, 257, 193, 147, 526, 173, 244, 286, 252, 353, 275, 442…
## $ E_PCI <int> -999, 29372, 22656, 20430, 24706, 22827, 27577, 24918, 2307…
## $ M_PCI <int> -999, 2306, 905, 1258, 758, 1482, 976, 1055, 854, 906, 925,…
## $ E_NOHSDP <int> 3669, 4204, 7861, 2141, 12620, 4383, 4833, 6016, 4530, 1040…
## $ M_NOHSDP <int> 426, 475, 727, 268, 766, 356, 429, 555, 426, 776, 464, 675,…
## $ E_AGE65 <int> 7083, 8050, 10233, 3806, 19386, 6409, 8375, 10502, 7650, 14…
## $ M_AGE65 <int> 25, 75, 91, 21, 119, 85, 75, 73, 100, 95, 85, 111, 99, 36, …
## $ E_AGE17 <int> 9318, 13369, 13468, 4566, 25196, 7006, 12161, 11608, 8190, …
## $ M_AGE17 <int> 14, 32, 53, 88, 66, 94, 0, 38, 85, 140, 39, 0, 92, 0, 0, 18…
## $ E_DISABL <int> 6280, 10465, 8114, 3492, 23598, 5570, 8788, 10086, 7703, 14…
## $ M_DISABL <int> 495, 729, 592, 370, 1086, 422, 488, 520, 416, 958, 664, 860…
## $ E_SNGPNT <int> 1330, 1586, 1437, 704, 4701, 1307, 1928, 1791, 1358, 2057, …
## $ M_SNGPNT <dbl> 285.0, 319.9, 267.2, 143.9, 464.0, 244.6, 272.8, 278.2, 195…
## $ E_MINRTY <int> 34397, 13788, 7413, 9641, 31675, 14954, 15065, 11630, 6205,…
## $ M_MINRTY <dbl> 145, 59, 229, 22, 34, 2, 56, 35, 25, 97, 187, 210, 87, 19, …
## $ E_LIMENG <int> 755, 426, 934, 93, 1076, 36, 664, 175, 129, 447, 461, 509, …
## $ M_LIMENG <dbl> 209.5, 205.9, 239.3, 137.4, 250.2, 87.9, 215.1, 127.6, 99.3…
## $ E_MUNIT <int> 67, 886, 211, 134, 1990, 679, 215, 587, 249, 959, 277, 895,…
## $ M_MUNIT <dbl> 37.1, 308.7, 104.2, 47.4, 303.0, 174.3, 99.9, 147.8, 86.8, …
## $ E_MOBILE <int> 7770, 4279, 6108, 2625, 7904, 2378, 3140, 2455, 4081, 8615,…
## $ M_MOBILE <int> 431, 469, 476, 212, 546, 261, 201, 310, 362, 672, 352, 493,…
## $ E_CROWD <int> 264, 299, 339, 119, 772, 404, 272, 83, 320, 555, 313, 422, …
## $ M_CROWD <dbl> 77.1, 142.3, 130.7, 57.7, 206.2, 148.9, 86.0, 54.7, 129.1, …
## $ E_NOVEH <int> 763, 1191, 856, 520, 2599, 989, 1189, 1322, 785, 1279, 1023…
## $ M_NOVEH <int> 160, 272, 201, 102, 331, 160, 218, 244, 170, 203, 187, 276,…
## $ E_GROUPQ <int> 654, 546, 543, 322, 3112, 512, 601, 436, 666, 1200, 1263, 4…
## $ M_GROUPQ <int> 142, 161, 117, 88, 436, 136, 81, 83, 175, 218, 214, 439, 16…
## $ EP_POV <dbl> -999.0, 15.4, 14.4, 23.5, 18.6, 16.6, 15.1, 16.3, 18.3, 15.…
## $ MP_POV <dbl> -999.0, 2.1, 1.7, 2.6, 1.2, 2.1, 1.4, 1.6, 1.9, 1.6, 1.8, 1…
## $ EP_UNEMP <dbl> -999.0, 4.2, 4.1, 6.7, 8.8, 5.0, 5.9, 5.9, 8.7, 5.5, 8.8, 6…
## $ MP_UNEMP <dbl> -999.0, 1.0, 0.8, 1.7, 1.0, 1.1, 1.0, 1.2, 1.5, 1.0, 1.4, 1…
## $ EP_PCI <dbl> -999, 29372, 22656, 20430, 24706, 22827, 27577, 24918, 2307…
## $ MP_PCI <dbl> -999, 2306, 905, 1258, 758, 1482, 976, 1055, 854, 906, 925,…
## $ EP_NOHSDP <dbl> 13.8, 11.3, 19.8, 15.4, 15.9, 18.6, 13.8, 15.6, 17.1, 18.2,…
## $ MP_NOHSDP <dbl> 1.6, 1.3, 1.8, 1.9, 1.0, 1.5, 1.2, 1.4, 1.6, 1.4, 1.4, 1.2,…
## $ EP_AGE65 <dbl> 18.0, 14.6, 17.8, 19.0, 16.8, 18.9, 16.3, 19.3, 20.5, 18.0,…
## $ MP_AGE65 <dbl> 0.1, 0.1, 0.2, 0.1, 0.1, 0.3, 0.1, 0.1, 0.3, 0.1, 0.2, 0.1,…
## $ EP_AGE17 <dbl> 23.7, 24.2, 23.4, 22.8, 21.9, 20.7, 23.7, 21.3, 21.9, 22.6,…
## $ MP_AGE17 <dbl> 0.0, 0.1, 0.1, 0.4, 0.1, 0.3, 0.0, 0.1, 0.2, 0.2, 0.1, 0.0,…
## $ EP_DISABL <dbl> 16.1, 19.3, 14.2, 17.7, 20.8, 16.7, 17.7, 18.7, 20.9, 17.6,…
## $ MP_DISABL <dbl> 1.3, 1.3, 1.0, 1.9, 0.9, 1.2, 1.0, 1.0, 1.1, 1.2, 1.4, 1.1,…
## $ EP_SNGPNT <dbl> 10.7, 7.5, 7.0, 10.5, 10.4, 9.7, 9.7, 8.2, 9.0, 6.7, 8.2, 9…
## $ MP_SNGPNT <dbl> 2.3, 1.5, 1.3, 2.1, 1.0, 1.8, 1.4, 1.3, 1.3, 0.8, 1.3, 1.4,…
## $ EP_MINRTY <dbl> 87.5, 25.0, 12.9, 48.1, 27.5, 44.2, 29.4, 21.3, 16.6, 7.9, …
## $ MP_MINRTY <dbl> 0.4, 0.1, 0.4, 0.1, 0.0, 0.0, 0.1, 0.1, 0.1, 0.1, 0.4, 0.3,…
## $ EP_LIMENG <dbl> 2.1, 0.8, 1.7, 0.5, 1.0, 0.1, 1.4, 0.3, 0.4, 0.6, 1.0, 0.7,…
## $ MP_LIMENG <dbl> 0.6, 0.4, 0.4, 0.7, 0.2, 0.3, 0.4, 0.2, 0.3, 0.2, 0.5, 0.3,…
## $ EP_MUNIT <dbl> 0.3, 3.8, 0.9, 1.3, 3.7, 4.0, 0.9, 2.2, 1.3, 2.5, 1.2, 2.6,…
## $ MP_MUNIT <dbl> 0.2, 1.3, 0.4, 0.5, 0.6, 1.0, 0.4, 0.6, 0.5, 0.5, 0.4, 0.8,…
## $ EP_MOBILE <dbl> 38.8, 18.4, 25.2, 26.2, 14.7, 14.0, 13.6, 9.3, 21.6, 22.9, …
## $ MP_MOBILE <dbl> 2.2, 2.0, 2.0, 2.1, 1.0, 1.5, 0.9, 1.2, 1.9, 1.8, 1.5, 1.5,…
## $ EP_CROWD <dbl> 2.1, 1.4, 1.6, 1.8, 1.7, 3.0, 1.4, 0.4, 2.1, 1.8, 1.7, 1.4,…
## $ MP_CROWD <dbl> 0.6, 0.7, 0.6, 0.9, 0.5, 1.1, 0.4, 0.3, 0.9, 0.5, 0.6, 0.6,…
## $ EP_NOVEH <dbl> 6.2, 5.6, 4.2, 7.8, 5.8, 7.3, 6.0, 6.1, 5.2, 4.2, 5.5, 4.7,…
## $ MP_NOVEH <dbl> 1.3, 1.3, 1.0, 1.5, 0.7, 1.2, 1.1, 1.1, 1.1, 0.7, 1.0, 0.9,…
## $ EP_GROUPQ <dbl> 1.7, 1.0, 0.9, 1.6, 2.7, 1.5, 1.2, 0.8, 1.8, 1.5, 2.6, 5.4,…
## $ MP_GROUPQ <dbl> 0.4, 0.3, 0.2, 0.4, 0.4, 0.4, 0.2, 0.2, 0.5, 0.3, 0.4, 0.5,…
## $ EPL_POV <dbl> -999.0000, 0.5401, 0.4723, 0.8860, 0.7322, 0.6258, 0.5229, …
## $ EPL_UNEMP <dbl> -999.0000, 0.2745, 0.2611, 0.6968, 0.8850, 0.4166, 0.5669, …
## $ EPL_PCI <dbl> -999.0000, 0.2860, 0.7561, 0.8879, 0.6076, 0.7465, 0.4022, …
## $ EPL_NOHSDP <dbl> 0.5922, 0.4397, 0.8405, 0.6753, 0.6988, 0.7966, 0.5922, 0.6…
## $ SPL_THEME1 <dbl> -999.0000, 1.5403, 2.3300, 3.1460, 2.9236, 2.5855, 2.0842, …
## $ RPL_THEME1 <dbl> -999.0000, 0.3631, 0.6143, 0.8455, 0.7866, 0.6901, 0.5373, …
## $ EPL_AGE65 <dbl> 0.4893, 0.1850, 0.4715, 0.5928, 0.3664, 0.5817, 0.3200, 0.6…
## $ EPL_AGE17 <dbl> 0.6826, 0.7529, 0.6406, 0.5578, 0.4323, 0.2846, 0.6826, 0.3…
## $ EPL_DISABL <dbl> 0.5610, 0.7905, 0.3763, 0.6845, 0.8564, 0.6074, 0.6845, 0.7…
## $ EPL_SNGPNT <dbl> 0.8383, 0.3792, 0.2961, 0.8185, 0.8077, 0.7316, 0.7316, 0.5…
## $ SPL_THEME2 <dbl> 2.5712, 2.1076, 1.7845, 2.6536, 2.4628, 2.2053, 2.4187, 2.2…
## $ RPL_THEME2 <dbl> 0.8758, 0.5810, 0.3187, 0.9077, 0.8303, 0.6609, 0.8067, 0.6…
## $ EPL_MINRTY <dbl> 0.9917, 0.6336, 0.4206, 0.8711, 0.6616, 0.8437, 0.6883, 0.5…
## $ EPL_LIMENG <dbl> 0.7740, 0.5113, 0.7170, 0.3582, 0.5791, 0.0707, 0.6718, 0.2…
## $ SPL_THEME3 <dbl> 1.7657, 1.1449, 1.1376, 1.2293, 1.2407, 0.9144, 1.3601, 0.8…
## $ RPL_THEME3 <dbl> 0.9268, 0.5947, 0.5915, 0.6447, 0.6507, 0.4597, 0.7256, 0.3…
## $ EPL_MUNIT <dbl> 0.0551, 0.6017, 0.1512, 0.2416, 0.5931, 0.6240, 0.1512, 0.4…
## $ EPL_MOBILE <dbl> 0.9869, 0.7408, 0.8816, 0.8940, 0.6396, 0.6192, 0.6065, 0.4…
## $ EPL_CROWD <dbl> 0.5498, 0.2964, 0.3703, 0.4457, 0.4056, 0.7622, 0.2964, 0.0…
## $ EPL_NOVEH <dbl> 0.5788, 0.4846, 0.2420, 0.7685, 0.5142, 0.7224, 0.5441, 0.5…
## $ EPL_GROUPQ <dbl> 0.4126, 0.1525, 0.1165, 0.3792, 0.6250, 0.3457, 0.2184, 0.0…
## $ SPL_THEME4 <dbl> 2.5832, 2.2760, 1.7616, 2.7290, 2.7775, 3.0735, 1.8166, 1.5…
## $ RPL_THEME4 <dbl> 0.5409, 0.3741, 0.1741, 0.6259, 0.6492, 0.7943, 0.1923, 0.1…
## $ SPL_THEMES <dbl> -999.0000, 7.0688, 7.0137, 9.7579, 9.4046, 8.7787, 7.6796, …
## $ RPL_THEMES <dbl> -999.0000, 0.4354, 0.4242, 0.8653, 0.8252, 0.7382, 0.5408, …
## $ F_POV <int> -999, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_UNEMP <int> -999, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_PCI <int> -999, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_NOHSDP <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_THEME1 <int> -999, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_AGE65 <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_AGE17 <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_DISABL <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_SNGPNT <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_THEME2 <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_MINRTY <int> 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_LIMENG <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_THEME3 <int> 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_MUNIT <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_MOBILE <int> 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_CROWD <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_NOVEH <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_GROUPQ <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_THEME4 <int> 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ F_TOTAL <int> -999, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ E_UNINSUR <int> 4160, 3875, 6303, 2005, 10686, 3613, 5004, 4353, 4116, 9631…
## $ M_UNINSUR <int> 588, 508, 732, 340, 796, 476, 540, 552, 506, 933, 600, 816,…
## $ EP_UNINSUR <dbl> 10.6, 7.1, 11.0, 10.2, 9.4, 10.8, 10.1, 8.1, 11.2, 11.8, 11…
## $ MP_UNINSUR <dbl> 1.5, 0.9, 1.3, 1.7, 0.7, 1.4, 1.1, 1.0, 1.4, 1.1, 1.3, 1.1,…
## $ E_DAYPOP <int> 32290, 37301, 40036, 17280, 117894, 27176, 44908, 56580, 35…dim(svi_2018)## [1] 3142 123The USGS glyphosate data has 15314 rows and 6 columns. The US County shapefile data has 3221 rows and 10 columns. The PRISM county-level data for 2013 has 1134420 rows and 6 columns. These data are connected by FIPS code, which is comprised of a two-digit state code and then a three digit county code. The USGS pesticide estimate data contains a low and high estimate (EPEST_LOW_KG and EPEST_HIGH_KG). The CDC SVI data has 3142 rows and 133 columns.
For the SVI, the four summary theme ranking variables, detailed in the Data Dictionary below, are:
Socioeconomic theme – RPL_THEME1 Household Composition and Disability – RPL_THEME2 Minority Status & Language – RPL_THEME3 Housing & Transportation – RPL_THEME4 The overall tract summary ranking variable is RPL_THEMES.
Variables that start with E are the estimates and M are the margins of error. Full documentation for SVI data can be found here: https://www.atsdr.cdc.gov/placeandhealth/svi/documentation/SVI_documentation_2014.html.
Data Cleaning/Wrangling
Create combined FIPS code (two-digit state plus three-digit county in format 01001).
usgs_gly_2013 <- usgs_gly_2013 %>% mutate(STATE_FIPS_CODE =
str_pad(as.character(usgs_gly_2013$STATE_FIPS_CODE),
2, pad = "0"),
COUNTY_FIPS_CODE =
str_pad(as.character(usgs_gly_2013$COUNTY_FIPS_CODE),
3, pad = "0"),
fips = as.factor(paste0(STATE_FIPS_CODE, COUNTY_FIPS_CODE)))
usgs_gly_2014 <- usgs_gly_2014 %>% mutate(STATE_FIPS_CODE =
str_pad(as.character(usgs_gly_2014$STATE_FIPS_CODE),
2, pad = "0"),
COUNTY_FIPS_CODE =
str_pad(as.character(usgs_gly_2014$COUNTY_FIPS_CODE),
3, pad = "0"),
fips = as.factor(paste0(STATE_FIPS_CODE, COUNTY_FIPS_CODE)))
usgs_gly_2015 <- usgs_gly_2015 %>% mutate(STATE_FIPS_CODE =
str_pad(as.character(usgs_gly_2015$STATE_FIPS_CODE),
2, pad = "0"),
COUNTY_FIPS_CODE =
str_pad(as.character(usgs_gly_2015$COUNTY_FIPS_CODE),
3, pad = "0"),
fips = as.factor(paste0(STATE_FIPS_CODE, COUNTY_FIPS_CODE)))
usgs_gly_2016 <- usgs_gly_2016 %>% mutate(STATE_FIPS_CODE =
str_pad(as.character(usgs_gly_2016$STATE_FIPS_CODE),
2, pad = "0"),
COUNTY_FIPS_CODE =
str_pad(as.character(usgs_gly_2016$COUNTY_FIPS_CODE),
3, pad = "0"),
fips = as.factor(paste0(STATE_FIPS_CODE, COUNTY_FIPS_CODE)))
usgs_gly_2017 <- usgs_gly_2017 %>% mutate(STATE_FIPS_CODE =
str_pad(as.character(usgs_gly_2017$STATE_FIPS_CODE),
2, pad = "0"),
COUNTY_FIPS_CODE =
str_pad(as.character(usgs_gly_2017$COUNTY_FIPS_CODE),
3, pad = "0"),
fips = as.factor(paste0(STATE_FIPS_CODE, COUNTY_FIPS_CODE)))
us_county_2013 <- us_county_2013 %>% mutate(fips =
as.factor(paste0(STATEFP, COUNTYFP)))
svi_2014 <- svi_2014 %>% mutate(fips =
str_pad(as.character(svi_2014$FIPS),
5, pad = "0"))
svi_2016 <- svi_2016 %>% mutate(fips =
str_pad(as.character(svi_2016$FIPS),
5, pad = "0"))
svi_2018 <- svi_2018 %>% mutate(fips =
str_pad(as.character(svi_2018$FIPS),
5, pad = "0"))#2013
prism_county_2013 <- prism_county_2013 %>% mutate(date=as.Date(date,
format="%d/%m/%Y")) %>%
select(-c(day,month, year))
group_by(prism_county_2013, fips) %>% tally()## # A tibble: 3,108 × 2
## fips n
## <fct> <int>
## 1 01001 365
## 2 01003 365
## 3 01005 365
## 4 01007 365
## 5 01009 365
## 6 01011 365
## 7 01013 365
## 8 01015 365
## 9 01017 365
## 10 01019 365
## # ℹ 3,098 more rows#got 3108 fips codes
prism_01001 <- prism_county_2013 %>% filter(fips=="01001")
quantile(prism_01001$tmean, probs = 0.85)## 85%
## 26.03194prism_12033 <- prism_county_2013 %>% filter(fips=="12033")
quantile(prism_12033$tmean, probs = 0.85)## 85%
## 26.83704prism_county_2013_85s <- prism_county_2013 %>% group_by(fips) %>% mutate(tmean85_fullyear = quantile(tmean, probs=0.85, na.rm=TRUE)) %>%
filter(between(date, as.Date('2013-04-30'), as.Date('2013-09-01'))) %>%
mutate(tmean85_spray = quantile(tmean, probs=0.85, na.rm=TRUE)) %>%
summarise(tmean85_fullyear=mean(tmean85_fullyear),
tmean85_spray=mean(tmean85_spray))
#2014
prism_county_2014 <- prism_county_2014 %>% mutate(date=as.Date(date,
format="%d/%m/%Y")) %>%
select(-c(day,month, year))
group_by(prism_county_2014, fips) %>% tally()## # A tibble: 3,108 × 2
## fips n
## <fct> <int>
## 1 01001 365
## 2 01003 365
## 3 01005 365
## 4 01007 365
## 5 01009 365
## 6 01011 365
## 7 01013 365
## 8 01015 365
## 9 01017 365
## 10 01019 365
## # ℹ 3,098 more rows#got 3108 fips codes
prism_01001 <- prism_county_2014 %>% filter(fips=="01001")
quantile(prism_01001$tmean, probs = 0.85)## 85%
## 26.47876prism_12033 <- prism_county_2014 %>% filter(fips=="12033")
quantile(prism_12033$tmean, probs = 0.85)## 85%
## 27.25865prism_county_2014_85s <- prism_county_2014 %>% group_by(fips) %>% mutate(tmean85_fullyear = quantile(tmean, probs=0.85, na.rm=TRUE)) %>%
filter(between(date, as.Date('2014-04-30'), as.Date('2014-09-01'))) %>%
mutate(tmean85_spray = quantile(tmean, probs=0.85, na.rm=TRUE)) %>%
summarise(tmean85_fullyear=mean(tmean85_fullyear),
tmean85_spray=mean(tmean85_spray))
#2015
prism_county_2015 <- prism_county_2015 %>% mutate(date=as.Date(date,
format="%d/%m/%Y")) %>%
select(-c(day,month, year))
group_by(prism_county_2015, fips) %>% tally()## # A tibble: 3,108 × 2
## fips n
## <fct> <int>
## 1 01001 365
## 2 01003 365
## 3 01005 365
## 4 01007 365
## 5 01009 365
## 6 01011 365
## 7 01013 365
## 8 01015 365
## 9 01017 365
## 10 01019 365
## # ℹ 3,098 more rows#got 3108 fips codes
prism_county_2015_85s <- prism_county_2015 %>% group_by(fips) %>% mutate(tmean85_fullyear = quantile(tmean, probs=0.85, na.rm=TRUE)) %>%
filter(between(date, as.Date('2015-04-30'), as.Date('2015-09-01'))) %>%
mutate(tmean85_spray = quantile(tmean, probs=0.85, na.rm=TRUE)) %>%
summarise(tmean85_fullyear=mean(tmean85_fullyear),
tmean85_spray=mean(tmean85_spray))
#2016
prism_county_2016 <- prism_county_2016 %>% mutate(date=as.Date(date,
format="%d/%m/%Y")) %>%
select(-c(day,month, year))
group_by(prism_county_2016, fips) %>% tally()## # A tibble: 3,108 × 2
## fips n
## <fct> <int>
## 1 01001 366
## 2 01003 366
## 3 01005 366
## 4 01007 366
## 5 01009 366
## 6 01011 366
## 7 01013 366
## 8 01015 366
## 9 01017 366
## 10 01019 366
## # ℹ 3,098 more rows#got 3108 fips codes
prism_county_2016_85s <- prism_county_2016 %>% group_by(fips) %>% mutate(tmean85_fullyear = quantile(tmean, probs=0.85, na.rm=TRUE)) %>%
filter(between(date, as.Date('2016-04-30'), as.Date('2016-09-01'))) %>%
mutate(tmean85_spray = quantile(tmean, probs=0.85, na.rm=TRUE)) %>%
summarise(tmean85_fullyear=mean(tmean85_fullyear),
tmean85_spray=mean(tmean85_spray))
#2017
prism_county_2017 <- prism_county_2017 %>% mutate(date=as.Date(date,
format="%d/%m/%Y")) %>%
select(-c(day,month, year))
group_by(prism_county_2017, fips) %>% tally()## # A tibble: 3,108 × 2
## fips n
## <fct> <int>
## 1 01001 365
## 2 01003 365
## 3 01005 365
## 4 01007 365
## 5 01009 365
## 6 01011 365
## 7 01013 365
## 8 01015 365
## 9 01017 365
## 10 01019 365
## # ℹ 3,098 more rows#got 3108 fips codes
prism_county_2017_85s <- prism_county_2017 %>% group_by(fips) %>% mutate(tmean85_fullyear = quantile(tmean, probs=0.85, na.rm=TRUE)) %>%
filter(between(date, as.Date('2017-04-30'), as.Date('2017-09-01'))) %>%
mutate(tmean85_spray = quantile(tmean, probs=0.85, na.rm=TRUE)) %>%
summarise(tmean85_fullyear=mean(tmean85_fullyear),
tmean85_spray=mean(tmean85_spray))svi_2014 <- svi_2014 %>% select(c(AFFGEOID, ST, STATE, ST_ABBR, COUNTY,fips,
RPL_THEME1, RPL_THEME2, RPL_THEME3,
RPL_THEME4, RPL_THEMES, Shape,
"Shape.STArea()","Shape.STLength()")) %>%
rename(SVI_socioecon=RPL_THEME1, SVI_housecompdis=RPL_THEME2,
SVI_minstatlang=RPL_THEME3, SVI_housing=RPL_THEME4,
SVI_overall=RPL_THEMES)
svi_2016 <- svi_2016 %>% select(c(ST, STATE, ST_ABBR, COUNTY,fips,
RPL_THEME1, RPL_THEME2, RPL_THEME3,
RPL_THEME4, RPL_THEMES)) %>%
rename(SVI_socioecon=RPL_THEME1, SVI_housecompdis=RPL_THEME2,
SVI_minstatlang=RPL_THEME3, SVI_housing=RPL_THEME4,
SVI_overall=RPL_THEMES)
#found a -999 value in svi_2018, filtering out
svi_2018 <- svi_2018 %>% select(c(ST, STATE, ST_ABBR, COUNTY,fips,
RPL_THEME1, RPL_THEME2, RPL_THEME3,
RPL_THEME4, RPL_THEMES)) %>%
rename(SVI_socioecon=RPL_THEME1, SVI_housecompdis=RPL_THEME2,
SVI_minstatlang=RPL_THEME3, SVI_housing=RPL_THEME4,
SVI_overall=RPL_THEMES) %>% filter(SVI_overall!=-999)
svi_2014 %>% summarize(mean_svi=mean(SVI_overall), sd_svi=sd(SVI_overall),
max_svi=max(SVI_overall), min_svi=min(SVI_overall),
n=length(unique(fips)))## mean_svi sd_svi max_svi min_svi n
## 1 0.4999994 0.2888131 1 0 3142svi_2016 %>% summarize(mean_svi=mean(SVI_overall), sd_svi=sd(SVI_overall),
max_svi=max(SVI_overall), min_svi=min(SVI_overall),
n=length(unique(fips)))## mean_svi sd_svi max_svi min_svi n
## 1 0.4999993 0.288813 1 0 3142svi_2018 %>% summarize(mean_svi=mean(SVI_overall), sd_svi=sd(SVI_overall),
max_svi=max(SVI_overall), min_svi=min(SVI_overall),
n=length(unique(fips)))## mean_svi sd_svi max_svi min_svi n
## 1 0.499993 0.2888138 1 0 3141merged_2013 <- left_join(us_county_2013, prism_county_2013_85s) %>%
left_join(usgs_gly_2013) %>% filter(STATEFP != "02" & STATEFP != "15" &
STATEFP != "72") %>%
mutate(gly_per_land_low = EPEST_LOW_KG/(ALAND/1e+6),
gly_per_land_high = EPEST_HIGH_KG/(ALAND/1e+6),
gly_per_totalarea_low = EPEST_LOW_KG/((ALAND+AWATER)/1e+6),
gly_per_totalarea_high = EPEST_HIGH_KG/((ALAND+AWATER)/1e+6)) %>%
left_join(svi_2014)## Joining with `by = join_by(fips)`
## Joining with `by = join_by(fips)`
## Joining with `by = join_by(AFFGEOID, fips)`merged_2014 <- left_join(us_county_2013, prism_county_2014_85s) %>%
left_join(usgs_gly_2014) %>% filter(STATEFP != "02" & STATEFP != "15" &
STATEFP != "72") %>%
mutate(gly_per_land_low = EPEST_LOW_KG/(ALAND/1e+6),
gly_per_land_high = EPEST_HIGH_KG/(ALAND/1e+6),
gly_per_totalarea_low = EPEST_LOW_KG/((ALAND+AWATER)/1e+6),
gly_per_totalarea_high = EPEST_HIGH_KG/((ALAND+AWATER)/1e+6)) %>%
left_join(svi_2014)## Joining with `by = join_by(fips)`
## Joining with `by = join_by(fips)`
## Joining with `by = join_by(AFFGEOID, fips)`merged_2015 <- left_join(us_county_2013, prism_county_2015_85s) %>%
left_join(usgs_gly_2015) %>% filter(STATEFP != "02" & STATEFP != "15" &
STATEFP != "72") %>%
mutate(gly_per_land_low = EPEST_LOW_KG/(ALAND/1e+6),
gly_per_land_high = EPEST_HIGH_KG/(ALAND/1e+6),
gly_per_totalarea_low = EPEST_LOW_KG/((ALAND+AWATER)/1e+6),
gly_per_totalarea_high = EPEST_HIGH_KG/((ALAND+AWATER)/1e+6)) %>%
left_join(svi_2016)## Joining with `by = join_by(fips)`
## Joining with `by = join_by(fips)`
## Joining with `by = join_by(fips)`merged_2016 <- left_join(us_county_2013, prism_county_2016_85s) %>%
left_join(usgs_gly_2016) %>% filter(STATEFP != "02" & STATEFP != "15" &
STATEFP != "72") %>%
mutate(gly_per_land_low = EPEST_LOW_KG/(ALAND/1e+6),
gly_per_land_high = EPEST_HIGH_KG/(ALAND/1e+6),
gly_per_totalarea_low = EPEST_LOW_KG/((ALAND+AWATER)/1e+6),
gly_per_totalarea_high = EPEST_HIGH_KG/((ALAND+AWATER)/1e+6)) %>%
left_join(svi_2016)## Joining with `by = join_by(fips)`
## Joining with `by = join_by(fips)`
## Joining with `by = join_by(fips)`merged_2017 <- left_join(us_county_2013, prism_county_2017_85s) %>%
left_join(usgs_gly_2017) %>% filter(STATEFP != "02" & STATEFP != "15" &
STATEFP != "72") %>%
mutate(gly_per_land_low = EPEST_LOW_KG/(ALAND/1e+6),
gly_per_land_high = EPEST_HIGH_KG/(ALAND/1e+6),
gly_per_totalarea_low = EPEST_LOW_KG/((ALAND+AWATER)/1e+6),
gly_per_totalarea_high = EPEST_HIGH_KG/((ALAND+AWATER)/1e+6)) %>%
left_join(svi_2018)## Joining with `by = join_by(fips)`
## Joining with `by = join_by(fips)`
## Joining with `by = join_by(fips)`#cleaned merged data so that it's just for conterminous U.S. (not AK which is 02, HI which is 15, or PR which is )compounddata_formap <- function(compound, year) {
usgs_data <- usgs_ests %>% filter(COMPOUND==compound) %>%
filter(YEAR==year) %>% mutate(STATE_FIPS_CODE =
str_pad(as.character(STATE_FIPS_CODE),
2, pad = "0"),
COUNTY_FIPS_CODE =
str_pad(as.character(COUNTY_FIPS_CODE),
3, pad = "0"),
fips = as.factor(paste0(STATE_FIPS_CODE,
COUNTY_FIPS_CODE)))
if (year==2013) {
us_county_data <- us_county_2013
prism_county_data <- prism_county_2013_85s
svi_data <- svi_2014
}
if (year==2014) {
prism_county_data <- prism_county_2014_85s
svi_data <- svi_2014
}
if (year==2015) {
prism_county_data <- prism_county_2015_85s
svi_data <- svi_2016
}
if (year==2016) {
prism_county_data <- prism_county_2016_85s
svi_data <- svi_2016
}
if (year==2017) {
prism_county_data <- prism_county_2017_85s
svi_data <- svi_2018
}
merged_data <- left_join(us_county_2013, prism_county_data) %>%
left_join(usgs_data) %>% filter(STATEFP != "02" & STATEFP != "15" &
STATEFP != "72") %>%
mutate(epest_per_land_low = EPEST_LOW_KG/(ALAND/1e+6),
epest_per_land_high = EPEST_HIGH_KG/(ALAND/1e+6),
epest_per_totalarea_low = EPEST_LOW_KG/((ALAND+AWATER)/1e+6),
epest_per_totalarea_high = EPEST_HIGH_KG/((ALAND+AWATER)/1e+6)) %>%
left_join(svi_data)
merged_data
}
hot_compounddata_formap <- function(compound, year, temp_c) {
compounddata_formap(compound, year) %>% filter(tmean85_spray>=temp_c)
}Data Visualization
Glyphosate Usage 2013-2017
plot_gly <- function(year, est_type) {
if (year==2013) {
data = merged_2013
}
if (year==2014) {
data = merged_2014
}
if (year==2015) {
data = merged_2015
}
if (year==2016) {
data = merged_2016
}
if (year==2017) {
data = merged_2017
}
if (est_type=="low") {
map <- ggplot(data) +
geom_sf(aes(fill = gly_per_land_low), linetype=0) +
ggtitle("Glyphosate per Land Area (Low Est.)", subtitle=year) +
theme_void(base_size = 14) +
theme(
plot.title = element_text(hjust = 0.5),
axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.ticks = element_blank(),
rect = element_blank()
) +
scale_fill_gradientn("Glyphosate usage \n(kg/km^2)",
colors = met.brewer("Tam"),
na.value = "grey50")
}
if (est_type=="high") {
map <- ggplot(data) +
geom_sf(aes(fill = gly_per_land_high), linetype=0) +
ggtitle("Glyphosate per Land Area (High Est.)", subtitle=year) +
theme_void(base_size = 14) +
theme(
plot.title = element_text(hjust = 0.5),
axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.ticks = element_blank(),
rect = element_blank()
) +
scale_fill_gradientn("Glyphosate usage \n(kg/km^2)",
colors = met.brewer("Tam"),
na.value = "grey50")
}
map
}
plot_gly(2013, "low")
plot_gly(2013, "high")
plot_gly(2014, "low")
plot_gly(2014, "high")
plot_gly(2015, "low")
plot_gly(2015, "high")
plot_gly(2016, "low")
plot_gly(2016, "high")
plot_gly(2017, "low")
plot_gly(2017, "high")
85th Percentile Mean Daily Temperatures 2013-2017
I created plots for temperature during spray season and during the full year.
2013
temp85_plot <- ggplot() +
geom_sf(data = merged_2013, aes(fill = tmean85_spray), linetype=0) +
ggtitle("85th Percentile Mean Temp. by County (2013)",
subtitle="Conterminous U.S.") + theme_void(base_size = 14) +
theme(
plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(colour="gray25", size=14,hjust = 0.5),
axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.ticks = element_blank(),
rect = element_blank()
) +
scale_fill_gradientn("Temperature (C)",
colors = rev(met.brewer("Hiroshige")),
na.value = "grey50")
temp85_plot_f_2013 <- ggplot() +
geom_sf(data = merged_2013, aes(fill = (tmean85_spray*9/5+32)), linetype=0) +
ggtitle("85th Percentile Mean Temp. by County (2013)",
subtitle="Conterminous U.S.") +
theme_void(base_size = 14) +
theme(
plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(colour="gray25", size=14,hjust = 0.5),
axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.ticks = element_blank(),
rect = element_blank()
) +
scale_fill_gradientn("Temperature (F)",
colors = rev(met.brewer("Hiroshige")),
na.value = "grey50")
plot(temp85_plot)
plot(temp85_plot_f_2013)
2014
temp85_plot <- ggplot() +
geom_sf(data = merged_2014, aes(fill = tmean85_spray), linetype=0) +
ggtitle("85th Percentile Mean Temp. by County (2014)",
subtitle="Conterminous U.S.") + theme_void(base_size = 14) +
theme(
plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(colour="gray25", size=14,hjust = 0.5),
axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.ticks = element_blank(),
rect = element_blank()
) +
scale_fill_gradientn("Temperature (C)",
colors = rev(met.brewer("Hiroshige")),
na.value = "grey50")
temp85_plot_f <- ggplot() +
geom_sf(data = merged_2014, aes(fill = (tmean85_spray*9/5+32)), linetype=0) +
ggtitle("85th Percentile Mean Temp. by County (2014)",
subtitle="Conterminous U.S.") +
theme_void(base_size = 14) +
theme(
plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(colour="gray25", size=14,hjust = 0.5),
axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.ticks = element_blank(),
rect = element_blank()
) +
scale_fill_gradientn("Temperature (F)",
colors = rev(met.brewer("Hiroshige")),
na.value = "grey50")
plot(temp85_plot)
plot(temp85_plot_f)
plot(temp85_plot)
plot_gly(2014, "high")
2015
temp85_plot <- ggplot() +
geom_sf(data = merged_2015, aes(fill = tmean85_spray), linetype=0) +
ggtitle("85th Percentile Mean Temp. by County (2015)",
subtitle="Conterminous U.S.") + theme_void(base_size = 14) +
theme(
plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(colour="gray25", size=14,hjust = 0.5),
axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.ticks = element_blank(),
rect = element_blank()
) +
scale_fill_gradientn("Temperature (C)",
colors = rev(met.brewer("Hiroshige")),
na.value = "grey50")
temp85_plot_f <- ggplot() +
geom_sf(data = merged_2015, aes(fill = (tmean85_spray*9/5+32)), linetype=0) +
ggtitle("85th Percentile Mean Temp. by County (2015)",
subtitle="Conterminous U.S.") +
theme_void(base_size = 14) +
theme(
plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(colour="gray25", size=14,hjust = 0.5),
axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.ticks = element_blank(),
rect = element_blank()
) +
scale_fill_gradientn("Temperature (F)",
colors = rev(met.brewer("Hiroshige")),
na.value = "grey50")
plot(temp85_plot)
plot(temp85_plot_f)
2016
temp85_plot <- ggplot() +
geom_sf(data = merged_2016, aes(fill = tmean85_spray), linetype=0) +
ggtitle("85th Percentile Mean Temp. by County (2016)",
subtitle="Conterminous U.S.") + theme_void(base_size = 14) +
theme(
plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(colour="gray25", size=14,hjust = 0.5),
axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.ticks = element_blank(),
rect = element_blank()
) +
scale_fill_gradientn("Temperature (C)",
colors = rev(met.brewer("Hiroshige")),
na.value = "grey50")
temp85_plot_f <- ggplot() +
geom_sf(data = merged_2016, aes(fill = (tmean85_spray*9/5+32)), linetype=0) +
ggtitle("85th Percentile Mean Temp. by County (2016)",
subtitle="Conterminous U.S.") +
theme_void(base_size = 14) +
theme(
plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(colour="gray25", size=14,hjust = 0.5),
axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.ticks = element_blank(),
rect = element_blank()
) +
scale_fill_gradientn("Temperature (F)",
colors = rev(met.brewer("Hiroshige")),
na.value = "grey50")
plot(temp85_plot)
plot(temp85_plot_f)
2017
temp85_plot <- ggplot() +
geom_sf(data = merged_2017, aes(fill = tmean85_spray), linetype=0) +
ggtitle("85th Percentile Mean Temp. by County (2017)",
subtitle="Conterminous U.S.") + theme_void(base_size = 14) +
theme(
plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(colour="gray25", size=14,hjust = 0.5),
axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.ticks = element_blank(),
rect = element_blank()
) +
scale_fill_gradientn("Temperature (C)",
colors = rev(met.brewer("Hiroshige")),
na.value = "grey50")
temp85_plot_f_2017 <- ggplot() +
geom_sf(data = merged_2017, aes(fill = (tmean85_spray*9/5+32)), linetype=0) +
ggtitle("85th Percentile Mean Temp. by County (2017)",
subtitle="Conterminous U.S.") +
theme_void(base_size = 14) +
theme(
plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(colour="gray25", size=14,hjust = 0.5),
axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.ticks = element_blank(),
rect = element_blank()
) +
scale_fill_gradientn("Temperature (F)",
colors = rev(met.brewer("Hiroshige")),
na.value = "grey50")
plot(temp85_plot)
plot(temp85_plot_f)
plot(temp85_plot)
plot_gly(2017, "high")
Higher SVI Counties
# Function to create SVI map
map_svi <- function(data, svi_var, year, percentile) {
# Filter data for the given year
filtered_data <- data %>%
filter(YEAR == year) %>%
filter(!!sym(svi_var) >= percentile) # Use rlang's !!sym to evaluate svi_var correctly
# Create the map
ggplot() +
# Layer for counties with high percentile (bold black outline)
geom_sf(
data = filtered_data,
aes(fill = !!sym(svi_var)), # Ensure the svi_var is evaluated correctly in aes
size = 2.0, # Thicker border
color = "black"
) +
# Layer for counties below the percentile (no fill or outline)
geom_sf(
data = data %>%
filter(YEAR == year) %>%
filter(!!sym(svi_var) < percentile),
aes(),
fill = NA,
color = NA
) +
scale_fill_viridis_c() + # Using viridis color scale for fill
labs(
title = paste("Map of", svi_var, "in", year),
fill = svi_var
) +
theme_minimal()
}
# Example of usage
# create_svi_map(data = my_data, svi_var = "SVI", year = 2023, percentile = 90, region_shapefile = "path/to/shapefile.shp")
#housing map for 2017
map_svi(merged_2017, "SVI_housing", 2017, 0.749)
#socioeconomic SVI 2017
map_svi(merged_2017, "SVI_socioecon", 2017, 0.749)
#household composition and disability 2017
map_svi(merged_2017,"SVI_housecompdis", 2017, 0.749)
#minority status and language 2017
map_svi(merged_2017, "SVI_minstatlang", 2017, 0.749)
#overall
map_svi(merged_2017, "SVI_overall", 2017, 0.749)
Data Analysis
I found which counties are hottest during the spray season (85th percentile of mean daily temperature 80 degrees F or hotter), and mapped the glyphosate usage for these counties.
Glyphosate Usage in Hot Counties 2013
all_leaflet_low <- function(compound, year) {
#get data
data <- compounddata_formap(compound, year)
map <- data %>%
st_transform(crs = "+init=epsg:4326") %>%
leaflet(width = "100%") %>%
addProviderTiles(provider = "CartoDB.Voyager") %>%
addPolygons(
stroke = FALSE,
smoothFactor = 0,
fillOpacity = 0.5,
fillColor= ~
colorNumeric("YlOrRd", epest_per_land_low)(epest_per_land_low)) %>%
addLegend(pal = colorNumeric("YlOrRd", data$epest_per_land_low), # Define color palette
values = data$epest_per_land_low, # Values to map
title = "Pesticide Usage (kg/km^2)", # Legend title
position = "bottomright", # Position of the legend
na.label = "No data"
)
map
}
all_leaflet_high <- function(compound, year) {
#get data
data <- compounddata_formap(compound, year)
map <- data %>%
st_transform(crs = "+init=epsg:4326") %>%
leaflet(width = "100%") %>%
addProviderTiles(provider = "CartoDB.Voyager") %>%
addPolygons(
stroke = FALSE,
smoothFactor = 0,
fillOpacity = 0.5,
fillColor= ~
colorNumeric("YlOrRd", epest_per_land_high)(epest_per_land_high)) %>%
addLegend(pal = colorNumeric("YlOrRd", data$epest_per_land_high), # Define color palette
values = data$epest_per_land_high, # Values to map
title = "Pesticide Usage (kg/km^2)", # Legend title
position = "bottomright", # Position of the legend
na.label = "No data"
)
map
}
hot_leaflet_low <- function(compound, year, temp_c) {
#get data
data <- hot_compounddata_formap(compound, year, temp_c)
map <- data %>%
st_transform(crs = "+init=epsg:4326") %>%
leaflet(width = "100%") %>%
addProviderTiles(provider = "CartoDB.Voyager") %>%
addPolygons(
stroke = FALSE,
smoothFactor = 0,
fillOpacity = 0.5,
fillColor= ~
colorNumeric("YlOrRd", epest_per_land_low)(epest_per_land_low)) %>%
addLegend(pal = colorNumeric("YlOrRd", data$epest_per_land_low), # Define color palette
values = data$epest_per_land_low, # Values to map
title = "Pesticide Usage (kg/km^2)", # Legend title
position = "bottomright", # Position of the legend
na.label = "No data"
)
map
}
hot_leaflet_high <- function(compound, year, temp_c) {
#get data
data <- hot_compounddata_formap(compound, year, temp_c)
map <- data %>%
st_transform(crs = "+init=epsg:4326") %>%
leaflet(width = "100%") %>%
addProviderTiles(provider = "CartoDB.Voyager") %>%
addPolygons(
stroke = FALSE,
smoothFactor = 0,
fillOpacity = 0.5,
fillColor= ~
colorNumeric("YlOrRd", epest_per_land_high)(epest_per_land_high)) %>%
addLegend(pal = colorNumeric("YlOrRd", data$epest_per_land_high), # Define color palette
values = data$epest_per_land_high, # Values to map
title = "Pesticide Usage (kg/km^2)", # Legend title
position = "bottomright", # Position of the legend
na.label = "No data"
)
map
}#find areas that have 85th percentile mean daily temp of 80 degrees F or higher
hot_merged_2013 <- merged_2013 %>% filter(tmean85_spray>=26.66)
hot_plot <- ggplot() +
geom_sf(data = hot_merged_2013, aes(fill = gly_per_land_low), linetype=0) +
ggtitle("Glyphosate per Land Area",
subtitle="Hot Counties") + theme_void(base_size = 14) +
theme(
plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(colour="gray25", size=14,hjust = 0.5),
axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.ticks = element_blank(),
rect = element_blank()
) +
scale_fill_gradientn("Glyphosate usage",
colors = met.brewer("Tam"),
na.value = "grey50")
plot(hot_plot)
hot_leaflet <- hot_merged_2013 %>%
st_transform(crs = "+init=epsg:4326") %>%
leaflet(width = "100%") %>%
addProviderTiles(provider = "CartoDB.Voyager") %>%
addPolygons(
stroke = FALSE,
smoothFactor = 0,
fillOpacity = 0.5,
fillColor= ~
colorNumeric("YlOrRd",gly_per_land_low)(gly_per_land_low)) %>%
addLegend(pal = colorNumeric("YlOrRd", hot_merged_2013$gly_per_land_low), # Define color palette
values = hot_merged_2013$gly_per_land_low, # Values to map
title = "Glyphosate Usage (kg/km^2)", # Legend title
position = "bottomright", # Position of the legend
na.label = "No data"
)## Warning in CPL_crs_from_input(x): GDAL Message 1: +init=epsg:XXXX syntax is
## deprecated. It might return a CRS with a non-EPSG compliant axis order.hot_leafletGlyphosate Usage in Hot Counties 2014
#find areas that have 85th percentile mean daily temp of 80 degrees F or higher
hot_merged_2014 <- merged_2014 %>% filter(tmean85_spray>=26.66)
hot_plot <- ggplot() +
geom_sf(data = hot_merged_2014, aes(fill = gly_per_land_low), linetype=0) +
ggtitle("Glyphosate per Land Area",
subtitle="Hot Counties") + theme_void(base_size = 14) +
theme(
plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(colour="gray25", size=14,hjust = 0.5),
axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.ticks = element_blank(),
rect = element_blank()
) +
scale_fill_gradientn("Glyphosate usage",
colors = met.brewer("Tam"),
na.value = "grey50")
plot(hot_plot)
hot_leaflet <- hot_merged_2014 %>%
st_transform(crs = "+init=epsg:4326") %>%
leaflet(width = "100%") %>%
addProviderTiles(provider = "CartoDB.Voyager") %>%
addPolygons(
stroke = FALSE,
smoothFactor = 0,
fillOpacity = 0.5,
fillColor= ~
colorNumeric("YlOrRd",gly_per_land_low)(gly_per_land_low)) %>%
addLegend(pal = colorNumeric("YlOrRd", hot_merged_2014$gly_per_land_low), # Define color palette
values = hot_merged_2014$gly_per_land_low, # Values to map
title = "Glyphosate Usage (kg/km^2)", # Legend title
position = "bottomright" # Position of the legend
)
hot_leafletGlyphosate Usage in Hot Counties 2015
#find areas that have 85th percentile mean daily temp of 80 degrees F or higher
hot_merged_2015 <- merged_2015 %>% filter(tmean85_spray>=26.66)
hot_plot <- ggplot() +
geom_sf(data = hot_merged_2015, aes(fill = gly_per_land_low), linetype=0) +
ggtitle("Glyphosate per Land Area",
subtitle="Hot Counties") + theme_void(base_size = 14) +
theme(
plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(colour="gray25", size=14,hjust = 0.5),
axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.ticks = element_blank(),
rect = element_blank()
) +
scale_fill_gradientn("Glyphosate usage",
colors = met.brewer("Tam"),
na.value = "grey50")
plot(hot_plot)
hot_leaflet <- hot_merged_2015 %>%
st_transform(crs = "+init=epsg:4326") %>%
leaflet(width = "100%") %>%
addProviderTiles(provider = "CartoDB.Voyager") %>%
addPolygons(
stroke = FALSE,
smoothFactor = 0,
fillOpacity = 0.5,
fillColor= ~
colorNumeric("YlOrRd",gly_per_land_low)(gly_per_land_low)) %>%
addLegend(pal = colorNumeric("YlOrRd", hot_merged_2015$gly_per_land_low), # Define color palette
values = hot_merged_2015$gly_per_land_low, # Values to map
title = "Glyphosate Usage (kg/km^2)", # Legend title
position = "bottomright" # Position of the legend
)
hot_leafletGlyphosate Usage in Hot Counties 2016
#find areas that have 85th percentile mean daily temp of 80 degrees F or higher
hot_merged_2016 <- merged_2016 %>% filter(tmean85_spray>=26.66)
hot_plot <- ggplot() +
geom_sf(data = hot_merged_2016, aes(fill = gly_per_land_low), linetype=0) +
ggtitle("Glyphosate per Land Area",
subtitle="Hot Counties") + theme_void(base_size = 14) +
theme(
plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(colour="gray25", size=14,hjust = 0.5),
axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.ticks = element_blank(),
rect = element_blank()
) +
scale_fill_gradientn("Glyphosate usage",
colors = met.brewer("Tam"),
na.value = "grey50")
plot(hot_plot)
hot_leaflet <- hot_merged_2016 %>%
st_transform(crs = "+init=epsg:4326") %>%
leaflet(width = "100%") %>%
addProviderTiles(provider = "CartoDB.Voyager") %>%
addPolygons(
stroke = FALSE,
smoothFactor = 0,
fillOpacity = 0.5,
fillColor= ~
colorNumeric("YlOrRd",gly_per_land_low)(gly_per_land_low)) %>%
addLegend(pal = colorNumeric("YlOrRd", hot_merged_2016$gly_per_land_low), # Define color palette
values = hot_merged_2016$gly_per_land_low, # Values to map
title = "Glyphosate Usage (kg/km^2)", # Legend title
position = "bottomright" # Position of the legend
)
hot_leafletGlyphosate Usage in Hot Counties 2017
#find areas that have 85th percentile mean daily temp of 80 degrees F or higher
hot_merged_2017 <- merged_2017 %>% filter(tmean85_spray>=26.66)
hot_plot <- ggplot() +
geom_sf(data = hot_merged_2017, aes(fill = gly_per_land_low), linetype=0) +
ggtitle("Glyphosate per Land Area",
subtitle="Hot Counties") + theme_void(base_size = 14) +
theme(
plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(colour="gray25", size=14,hjust = 0.5),
axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.ticks = element_blank(),
rect = element_blank()
) +
scale_fill_gradientn("Glyphosate usage",
colors = met.brewer("Tam"),
na.value = "grey50")
plot(hot_plot)
hot_leaflet <- hot_merged_2017 %>%
st_transform(crs = "+init=epsg:4326") %>%
leaflet(width = "100%") %>%
addProviderTiles(provider = "CartoDB.Voyager") %>%
addPolygons(
stroke = FALSE,
smoothFactor = 0,
fillOpacity = 0.5,
fillColor= ~
colorNumeric("YlOrRd",gly_per_land_low)(gly_per_land_low)) %>%
addLegend(pal = colorNumeric("YlOrRd", hot_merged_2017$gly_per_land_low), # Define color palette
values = hot_merged_2017$gly_per_land_low, # Values to map
title = "Glyphosate Usage (kg/km^2)", # Legend title
position = "bottomright" # Position of the legend
)
hot_leafletChange in Hot Counties Over Time
hot_merged_all <- list(hot_merged_2013, hot_merged_2014, hot_merged_2015,
hot_merged_2016, hot_merged_2017) %>%
rbindlist(fill=TRUE) %>% mutate(fips=as.numeric(fips))
hot_merged_sums <- hot_merged_all %>% group_by(YEAR) %>%
summarise(states=length(unique(STATEFP)),
counties=length(unique(fips)),
tmean85_fullyear = mean(tmean85_fullyear),
tmean85_spray = mean(tmean85_spray),
gly_low = sum(EPEST_LOW_KG, na.rm=TRUE),
gly_high = sum(EPEST_HIGH_KG)) %>%
filter(!is.na(YEAR))
hot_gly_use <- ggplot(hot_merged_sums) +
aes(x = YEAR, y = gly_low) +
geom_area() +
geom_line(aes(x = YEAR, y = counties), colour = "#461711") +
labs(
x = "Year",
y = "Glyphosate usage (kg)",
title = "Glyphosate Usage 2013-2017",
subtitle = "Hot Counties"
) +
theme_bw() +
theme(
plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5)
)
plot(hot_gly_use)
Characterizing population in these areas
#put together simple summary table that counts how many people are in these hot counties. Include SVI & subscores.
hot_merged_2013 %>% group_by(STATEFP) %>%
# I don't actually have data on population!#Shiny App Good resources for helping learn Shiny: https://mastering-shiny.org/scaling-functions.html?q=map#functional-programming and example code creating an interactive map using Shiny: https://github.com/fverkroost/RStudio-Blogs/blob/master/interactive_worldmap_shiny_app.R These could also be useful: https://forum.posit.co/t/help-creating-a-simple-leaflet-interactive-plot-in-r-shiny/182374/2 https://rviews.rstudio.com/2019/10/09/building-interactive-world-maps-in-shiny/
# could create a binary data structure where there are 5 columns for each year and 448 rows for the compounds
#useful functions are compounddata_formap() and hot_compounddata_formap
ui <- fluidPage(
# Application title
titlePanel("Pesticides and Heat (2013-2017)"),
fluidRow(
column(3, selectizeInput("pesticide", label = "Pesticide", choices = usgs_compounds)),
column(3, selectInput("year", label = "Year", choices = c(2013, 2014, 2015, 2016, 2017))),
column(4, numericInput('temp', label="85th %ile temp (C) in spray season", value=25)),
column(3, checkboxInput('hot_counties', label= "Hottest counties only"))),
# Application layout
leafletOutput("map")
)
server = function(input, output, session) {
# render the map
output$map <- renderLeaflet({
if(input$hot_counties==TRUE) {
hot_leaflet_low(input$pesticide, input$year, input$temp)
}
else {
all_leaflet_low(input$pesticide, input$year)
}
})
}
shinyApp(ui, server)
# look at Ag Census and ag worker density
# flexdash flexdashboard
# are there systematic patterns of more being applied when it is hot?
#make scatterplots of heat x pesticide where each location is an observationsaveRDS(usgs_ests, file="usgs_ests.rds")
saveRDS(usgs_compounds, file="usgs_compounds.rds")
saveRDS(us_county_2013, file="us_county_data.rds")
saveRDS(prism_county_2013_85s, file="prism_county_2013_85s.rds")
saveRDS(prism_county_2014_85s, file="prism_county_2014_85s.rds")
saveRDS(prism_county_2015_85s, file="prism_county_2015_85s.rds")
saveRDS(prism_county_2016_85s, file="prism_county_2016_85s.rds")
saveRDS(prism_county_2017_85s, file="prism_county_2017_85s.rds")
saveRDS(svi_2014, file="svi_2014.rds")
saveRDS(svi_2016, file="svi_2016.rds")
saveRDS(svi_2018, file="svi_2018.rds")Functions
#-----functions----
# Show the names of all functions defined in the .Rmd
# (e.g. loaded in the environment)
lsf.str()## all_leaflet_high : function (compound, year)
## all_leaflet_low : function (compound, year)
## compounddata_formap : function (compound, year)
## hot_compounddata_formap : function (compound, year, temp_c)
## hot_leaflet_high : function (compound, year, temp_c)
## hot_leaflet_low : function (compound, year, temp_c)
## map_svi : function (data, svi_var, year, percentile)
## plot_gly : function (year, est_type)# Show the definitions of all functions loaded into the current environment
lsf.str() %>% set_names() %>% map(get, .GlobalEnv)## $all_leaflet_high
## <srcref: file "" chars 26:21 to 49:1>
##
## $all_leaflet_low
## <srcref: file "" chars 1:20 to 24:1>
##
## $compounddata_formap
## <srcref: file "" chars 1:24 to 48:1>
##
## $hot_compounddata_formap
## <srcref: file "" chars 50:28 to 52:1>
##
## $hot_leaflet_high
## <srcref: file "" chars 76:21 to 99:1>
##
## $hot_leaflet_low
## <srcref: file "" chars 51:20 to 74:1>
##
## $map_svi
## <srcref: file "" chars 2:12 to 32:1>
## <bytecode: 0x55b3408941a8>
##
## $plot_gly
## <srcref: file "" chars 1:13 to 52:1>
## <bytecode: 0x55b342a96940>Versions
Information about the version of R and the packages using are below:
sessionInfo()## R version 4.3.3 (2024-02-29)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 20.04.6 LTS
##
## Matrix products: default
## BLAS: /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.9.0
## LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.9.0
##
## locale:
## [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
## [3] LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8
## [5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8
## [7] LC_PAPER=en_US.UTF-8 LC_NAME=C
## [9] LC_ADDRESS=C LC_TELEPHONE=C
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
##
## time zone: Etc/UTC
## tzcode source: system (glibc)
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] shiny_1.9.1 leaflet_2.2.2 terra_1.7-78 data.table_1.16.2
## [5] sf_1.0-18 here_1.0.1 ggpubr_0.6.0 MetBrewer_0.2.0
## [9] lubridate_1.9.3 forcats_1.0.0 stringr_1.5.1 dplyr_1.1.4
## [13] purrr_1.0.2 readr_2.1.5 tibble_3.2.1 ggplot2_3.5.1
## [17] tidyverse_2.0.0 tidyr_1.3.1 RColorBrewer_1.1-3
##
## loaded via a namespace (and not attached):
## [1] gtable_0.3.5 xfun_0.48 bslib_0.8.0
## [4] htmlwidgets_1.6.4 rstatix_0.7.2 tzdb_0.4.0
## [7] leaflet.providers_2.0.0 crosstalk_1.2.1 vctrs_0.6.5
## [10] tools_4.3.3 generics_0.1.3 proxy_0.4-27
## [13] fansi_1.0.6 highr_0.11 pacman_0.5.1
## [16] pkgconfig_2.0.3 KernSmooth_2.23-22 lifecycle_1.0.4
## [19] farver_2.1.2 compiler_4.3.3 munsell_0.5.1
## [22] codetools_0.2-20 carData_3.0-5 httpuv_1.6.15
## [25] htmltools_0.5.8.1 class_7.3-22 sass_0.4.9
## [28] yaml_2.3.10 Formula_1.2-5 later_1.3.2
## [31] pillar_1.9.0 car_3.1-3 jquerylib_0.1.4
## [34] classInt_0.4-10 cachem_1.1.0 abind_1.4-8
## [37] mime_0.12 tidyselect_1.2.1 digest_0.6.37
## [40] stringi_1.8.4 labeling_0.4.3 rprojroot_2.0.4
## [43] fastmap_1.2.0 grid_4.3.3 colorspace_2.1-1
## [46] cli_3.6.3 magrittr_2.0.3 utf8_1.2.4
## [49] broom_1.0.7 e1071_1.7-16 withr_3.0.1
## [52] promises_1.3.0 scales_1.3.0 backports_1.5.0
## [55] timechange_0.3.0 rmarkdown_2.28 ggsignif_0.6.4
## [58] hms_1.1.3 evaluate_1.0.1 knitr_1.48
## [61] viridisLite_0.4.2 rlang_1.1.4 Rcpp_1.0.13
## [64] xtable_1.8-4 glue_1.8.0 DBI_1.2.3
## [67] rstudioapi_0.16.0 jsonlite_1.8.9 R6_2.5.1
## [70] units_0.8-5